New OctQuad Method

cmd/benchmark/main.go

@@ -38,7 +38,7 @@ func (ds dataset) WriteCSV(out io.Writer) (int, error) {
 
 		runtimeOut := "NA"
 		if e.duration != nil {
-			runtimeOut = fmt.Sprintf("%s", *e.duration)
+			runtimeOut = fmt.Sprintf("%d", e.duration.Milliseconds())
 		}
 
 		avgErrOut := "NA"
@@ -140,6 +140,30 @@ func (coarse24w coarse24Writer) method() string                 { return "coarse
 func (coarse24w coarse24Writer) pack(v vector.Vector3) []byte   { return unitpacking.PackCoarse24(v) }
 func (coarse24w coarse24Writer) unpack(b []byte) vector.Vector3 { return unitpacking.UnpackCoarse24(b) }
 
+type octQuad16Writer struct{ out io.Writer }
+
+func (octQuad16w octQuad16Writer) method() string               { return "octquad16" }
+func (octQuad16w octQuad16Writer) pack(v vector.Vector3) []byte { return unitpacking.PackOctQuad16(v) }
+func (octQuad16w octQuad16Writer) unpack(b []byte) vector.Vector3 {
+	return unitpacking.UnpackOctQuad16(b)
+}
+
+type octQuad24Writer struct{ out io.Writer }
+
+func (octQuad24w octQuad24Writer) method() string               { return "octquad24" }
+func (octQuad24w octQuad24Writer) pack(v vector.Vector3) []byte { return unitpacking.PackOctQuad24(v) }
+func (octQuad24w octQuad24Writer) unpack(b []byte) vector.Vector3 {
+	return unitpacking.UnpackOctQuad24(b)
+}
+
+type octQuad32Writer struct{ out io.Writer }
+
+func (octQuad32w octQuad32Writer) method() string               { return "octquad32" }
+func (octQuad32w octQuad32Writer) pack(v vector.Vector3) []byte { return unitpacking.PackOctQuad32(v) }
+func (octQuad32w octQuad32Writer) unpack(b []byte) vector.Vector3 {
+	return unitpacking.UnpackOctQuad32(b)
+}
+
 type oct16Writer struct{ out io.Writer }
 
 func (oct16w oct16Writer) method() string                 { return "oct16" }
@@ -246,7 +270,6 @@ func getDatasetPathsFromDir(dir string) ([]string, error) {
 }
 
 func calcFlatNormals(m mango.Mesh) []vector.Vector3 {
-
 	normals := make([]vector.Vector3, len(m.Vertices()))
 	for i := range normals {
 		normals[i] = vector.Vector3One()
@@ -442,7 +465,6 @@ func runbaseline(unitVectors []vector.Vector3) runResultEntry {
 }
 
 func writeObj(mesh mango.Mesh, normals []vector.Vector3, out io.Writer) error {
-
 	for _, n := range mesh.Vertices() {
 		_, err := fmt.Fprintf(out, "v %f %f %f\n", n.X(), n.Y(), n.Z())
 		if err != nil {
@@ -493,6 +515,9 @@ func main() {
 		oct16Writer{os.Stdout},
 		oct24Writer{os.Stdout},
 		oct32Writer{os.Stdout},
+		octQuad16Writer{os.Stdout},
+		octQuad24Writer{os.Stdout},
+		octQuad32Writer{os.Stdout},
 	}
 
 	if writeCSV {
@@ -546,5 +571,4 @@ func main() {
 			panic(err)
 		}
 	}
-
 }

cmd/example/main.go

@@ -1,6 +1,10 @@
 package main
 
 import (
+	"compress/flate"
+	"image"
+	"image/color"
+	"image/png"
 	"math/rand"
 	"os"
 
@@ -9,9 +13,16 @@ import (
 )
 
 func main() {
-	numOfVectors := 10000000
+
+	width := 512
+	height := 512
+	numOfVectors := width * height
 	unitVectors := make([]vector.Vector3, numOfVectors)
 
+	upLeft := image.Point{0, 0}
+	lowRight := image.Point{width, height}
+	img := image.NewRGBA(image.Rectangle{upLeft, lowRight})
+
 	// Generate a bunch of unit vectors
 	for i := 0; i < numOfVectors; i++ {
 		unitVectors[i] = vector.NewVector3(
@@ -26,8 +37,29 @@ func main() {
 		panic(err)
 	}
 
+	comressedWriter, err := flate.NewWriter(out, 9)
+	if err != nil {
+		panic(err)
+	}
+
 	// Write out unit vectors in packed format
-	for _, unitVector := range unitVectors {
-		out.Write(unitpacking.PackOct24(unitVector))
+	for x := 0; x < width; x++ {
+		for y := 0; y < height; y++ {
+			data := unitpacking.PackOct24(unitVectors[y+(x*width)])
+			comressedWriter.Write(data)
+			img.Set(x, y, color.RGBA{
+				R: data[0],
+				G: data[1],
+				B: data[2],
+				A: 255,
+			})
+		}
 	}
+
+	comressedWriter.Flush()
+
+	// Encode as PNG.
+	f, _ := os.Create("image.png")
+	png.Encode(f, img)
+
 }

unitpacking/data_test.go

@@ -19,3 +19,16 @@ var testVectors []vector.Vector3 = []vector.Vector3{
 	vector.NewVector3(-0.997605826445425, 0.06365823804882093, -0.027023022974122023),
 	vector.NewVector3(0.7180684556508264, -0.6958747397502424, -0.011663600506254649),
 }
+
+var quadTestVectors []vector.Vector2 = []vector.Vector2{
+	vector.NewVector2(1, 1),
+	vector.NewVector2(-1, -1),
+	vector.NewVector2(1, -1),
+	vector.NewVector2(-1, 1),
+	vector.NewVector2(0, 0),
+	vector.NewVector2(0.5, 0.5),
+	vector.NewVector2(0.5, 0.25),
+	vector.NewVector2(0, 0.25),
+	vector.NewVector2(0, -0.25),
+	vector.NewVector2(0.12, -.94),
+}

unitpacking/oct.go

@@ -31,7 +31,7 @@ func multVect(a, b vector.Vector2) vector.Vector2 {
 // PackOct32 maps a unit vector to a 2D UV of a octahedron, and then writes the
 // 2D coordinates to 4 bytes, 2 bytes per coordinate.
 func PackOct32(v vector.Vector3) []byte {
-	uvCords := MapToOctUV(v)
+	uvCords := MapToOctUVPrecise(v, 32)
 
 	// 2 ^ 16 = 65,536;
 	x := uint(math.Floor(uvCords.X()*32767) + 32768)
@@ -47,22 +47,22 @@ func PackOct32(v vector.Vector3) []byte {
 }
 
 // UnpackOct32 reads in two 16bit numbers and converts from 2D octahedron UV to
-// 3D unit sphere coordinates
+// 3D unit sphere coordinates.
 func UnpackOct32(b []byte) vector.Vector3 {
 	everything := uint(b[0]) | (uint(b[1]) << 8) | (uint(b[2]) << 16) | (uint(b[3]) << 24)
 	rawY := (int)((everything) & 0xFFFF)
 	rawX := (int)(everything >> 16)
 
-	cleanedX := clamp((float64(rawX)-32768.0)/32767.0, -1.0, 1.0)
-	cleanedY := clamp((float64(rawY)-32768.0)/32767.0, -1.0, 1.0)
+	cleanedX := Clamp((float64(rawX)-32768.0)/32767.0, -1.0, 1.0)
+	cleanedY := Clamp((float64(rawY)-32768.0)/32767.0, -1.0, 1.0)
 
 	return FromOctUV(vector.NewVector2(cleanedX, cleanedY))
 }
 
 // PackOct24 maps a unit vector to a 2D UV of a octahedron, and then writes the
 // 2D coordinates to 3 bytes, 12bits per coordinate.
 func PackOct24(v vector.Vector3) []byte {
-	uvCords := MapToOctUV(v)
+	uvCords := MapToOctUVPrecise(v, 24)
 
 	// 2 ^ 12 = 4,096;
 	x := uint(math.Floor(uvCords.X()*2047) + 2048)
@@ -77,22 +77,22 @@ func PackOct24(v vector.Vector3) []byte {
 }
 
 // UnpackOct24 reads in two 12bit numbers and converts from 2D octahedron UV to
-// 3D unit sphere coordinates
+// 3D unit sphere coordinates.
 func UnpackOct24(b []byte) vector.Vector3 {
 	everything := uint(b[0]) | (uint(b[1]) << 8) | (uint(b[2]) << 16)
 	rawY := (int)((everything) & 0b111111111111)
 	rawX := (int)(everything >> 12)
 
-	cleanedX := clamp((float64(rawX)-2048.0)/2047.0, -1.0, 1.0)
-	cleanedY := clamp((float64(rawY)-2048.0)/2047.0, -1.0, 1.0)
+	cleanedX := Clamp((float64(rawX)-2048.0)/2047.0, -1.0, 1.0)
+	cleanedY := Clamp((float64(rawY)-2048.0)/2047.0, -1.0, 1.0)
 
 	return FromOctUV(vector.NewVector2(cleanedX, cleanedY))
 }
 
 // PackOct16 maps a unit vector to a 2D UV of a octahedron, and then writes the
 // 2D coordinates to 2 bytes, 8bits per coordinate.
 func PackOct16(v vector.Vector3) []byte {
-	uvCords := MapToOctUV(v)
+	uvCords := MapToOctUVPrecise(v, 16)
 
 	// 2 ^ 8 = 256;
 	x := uint(math.Floor(uvCords.X()*127) + 128)
@@ -106,23 +106,63 @@ func PackOct16(v vector.Vector3) []byte {
 }
 
 // UnpackOct16 reads in two 8bit numbers and converts from 2D octahedron UV to
-// 3D unit sphere coordinates
+// 3D unit sphere coordinates.
 func UnpackOct16(b []byte) vector.Vector3 {
 	everything := uint(b[0]) | (uint(b[1]) << 8)
 	rawY := (int)((everything) & 0b11111111)
 	rawX := (int)(everything >> 8)
 
-	cleanedX := clamp((float64(rawX)-128.0)/127.0, -1.0, 1.0)
-	cleanedY := clamp((float64(rawY)-128.0)/127.0, -1.0, 1.0)
+	cleanedX := Clamp((float64(rawX)-128.0)/127.0, -1.0, 1.0)
+	cleanedY := Clamp((float64(rawY)-128.0)/127.0, -1.0, 1.0)
 
 	return FromOctUV(vector.NewVector2(cleanedX, cleanedY))
 }
 
-// MapToOctUV converts a 3D sphere's coordinates to a 2D octahedron UV
+// MapToOctUVPrecise brute force finds an optimal UV coordinate that minimizes
+// rounding error.
+func MapToOctUVPrecise(v vector.Vector3, n int) vector.Vector2 {
+	s := MapToOctUV(v) // Remap to the square
+
+	// Each snorm’s max value interpreted as an integer,
+	// e.g., 127.0 for snorm8
+	M := float64(int(1)<<((n/2)-1)) - 1.0
+
+	// Remap components to snorm(n/2) precision...with floor instead
+	// of round (see equation 1)
+	s = floorVec2(clampVec2(s, -1.0, 1.0).MultByConstant(M)).MultByConstant(1.0 / M)
+	bestRepresentation := s
+	highestCosine := FromOctUV(s).Dot(v)
+
+	// Test all combinations of floor and ceil and keep the best.
+	// Note that at +/- 1, this will exit the square... but that
+	// will be a worse encoding and never win.
+	for i := 0; i <= 1; i++ {
+		for j := 0; j <= 1; j++ {
+			// This branch will be evaluated at compile time
+			if (i != 0) || (j != 0) {
+				// Offset the bit pattern (which is stored in floating
+				// point!) to effectively change the rounding mode
+				// (when i or j is 0: floor, when it is one: ceiling)
+				candidate := vector.NewVector2(float64(i), float64(j)).MultByConstant(1 / M).Add(s)
+				cosine := FromOctUV(candidate).Dot(v)
+				if cosine > highestCosine {
+					bestRepresentation = candidate
+					highestCosine = cosine
+				}
+			}
+		}
+	}
+
+	return bestRepresentation
+}
+
+// MapToOctUV converts a 3D sphere's coordinates to a 2D octahedron UV.
 func MapToOctUV(v vector.Vector3) vector.Vector2 {
 	// Project the sphere onto the octahedron, and then onto the xy plane
 	// vec2 p = v.xy * (1.0 / (abs(v.x) + abs(v.y) + abs(v.z)));
-	p := vector.NewVector2(v.X(), v.Y()).MultByConstant(1.0 / (math.Abs(v.X()) + math.Abs(v.Y()) + math.Abs(v.Z())))
+	p := vector.
+		NewVector2(v.X(), v.Y()).
+		MultByConstant(1.0 / (math.Abs(v.X()) + math.Abs(v.Y()) + math.Abs(v.Z())))
 	if v.Z() > 0 {
 		return p
 	}

unitpacking/oct_quad.go

@@ -0,0 +1,39 @@
+package unitpacking
+
+import "github.com/EliCDavis/vector"
+
+// PackOctQuad16 maps a unit vector to a 2D UV of a octahedron, and then
+// encodes the 2D coordinates inside a quad tree.
+func PackOctQuad16(v vector.Vector3) []byte {
+	return Vec2ToTwoByteQuad(MapToOctUV(v))
+}
+
+// UnpackOctQuad16 builds a 2D coordinate from the encoded quadtree and then
+// converts the 2D octahedron UV to 3D unit sphere coordinates.
+func UnpackOctQuad16(b []byte) vector.Vector3 {
+	return FromOctUV(TwoByteQuadToVec2(b))
+}
+
+// PackOctQuad24 maps a unit vector to a 2D UV of a octahedron, and then
+// encodes the 2D coordinates inside a quad tree.
+func PackOctQuad24(v vector.Vector3) []byte {
+	return Vec2ToThreeByteQuad(MapToOctUV(v))
+}
+
+// UnpackOctQuad24 builds a 2D coordinate from the encoded quadtree and then
+// converts the 2D octahedron UV to 3D unit sphere coordinates.
+func UnpackOctQuad24(b []byte) vector.Vector3 {
+	return FromOctUV(ThreeByteQuadToVec2(b))
+}
+
+// PackOctQuad32 maps a unit vector to a 2D UV of a octahedron, and then
+// encodes the 2D coordinates inside a quad tree.
+func PackOctQuad32(v vector.Vector3) []byte {
+	return Vec2ToFourByteQuad(MapToOctUV(v))
+}
+
+// UnpackOctQuad32 builds a 2D coordinate from the encoded quadtree and then
+// converts the 2D octahedron UV to 3D unit sphere coordinates.
+func UnpackOctQuad32(b []byte) vector.Vector3 {
+	return FromOctUV(FourByteQuadToVec2(b))
+}

unitpacking/oct_quad_test.go

@@ -0,0 +1,60 @@
+package unitpacking_test
+
+import (
+	"fmt"
+	"testing"
+
+	"github.com/recolude/unitpacking/unitpacking"
+	"github.com/stretchr/testify/assert"
+)
+
+func TestOctQuad16(t *testing.T) {
+	for _, tc := range testVectors {
+		unit := tc.Normalized()
+		name := fmt.Sprintf("%.2f,%.2f,%.2f", unit.X(), unit.Y(), unit.Z())
+
+		t.Run(name, func(t *testing.T) {
+			packed := unitpacking.PackOctQuad16(unit)
+			assert.Len(t, packed, 2)
+			unpacked := unitpacking.UnpackOctQuad16(packed)
+
+			assert.InDelta(t, unit.X(), unpacked.X(), 0.02, "X components not equal: %.2f != %.2f", unit.X(), unpacked.X())
+			assert.InDelta(t, unit.Y(), unpacked.Y(), 0.02, "Y components not equal: %.2f != %.2f", unit.Y(), unpacked.Y())
+			assert.InDelta(t, unit.Z(), unpacked.Z(), 0.02, "Z components not equal: %.2f != %.2f", unit.Z(), unpacked.Z())
+		})
+	}
+}
+
+func TestOctQuad24(t *testing.T) {
+	for _, tc := range testVectors {
+		unit := tc.Normalized()
+		name := fmt.Sprintf("%.2f,%.2f,%.2f", unit.X(), unit.Y(), unit.Z())
+
+		t.Run(name, func(t *testing.T) {
+			packed := unitpacking.PackOctQuad24(unit)
+			assert.Len(t, packed, 3)
+			unpacked := unitpacking.UnpackOctQuad24(packed)
+
+			assert.InDelta(t, unit.X(), unpacked.X(), 0.02, "X components not equal: %.2f != %.2f", unit.X(), unpacked.X())
+			assert.InDelta(t, unit.Y(), unpacked.Y(), 0.02, "Y components not equal: %.2f != %.2f", unit.Y(), unpacked.Y())
+			assert.InDelta(t, unit.Z(), unpacked.Z(), 0.02, "Z components not equal: %.2f != %.2f", unit.Z(), unpacked.Z())
+		})
+	}
+}
+
+func TestOctQuad32(t *testing.T) {
+	for _, tc := range testVectors {
+		unit := tc.Normalized()
+		name := fmt.Sprintf("%.2f,%.2f,%.2f", unit.X(), unit.Y(), unit.Z())
+
+		t.Run(name, func(t *testing.T) {
+			packed := unitpacking.PackOctQuad32(unit)
+			assert.Len(t, packed, 4)
+			unpacked := unitpacking.UnpackOctQuad32(packed)
+
+			assert.InDelta(t, unit.X(), unpacked.X(), 0.00005, "X components not equal: %.2f != %.2f", unit.X(), unpacked.X())
+			assert.InDelta(t, unit.Y(), unpacked.Y(), 0.00005, "Y components not equal: %.2f != %.2f", unit.Y(), unpacked.Y())
+			assert.InDelta(t, unit.Z(), unpacked.Z(), 0.00005, "Z components not equal: %.2f != %.2f", unit.Z(), unpacked.Z())
+		})
+	}
+}