From ccd3337b79d871c967e7145b6396b96ffa3923ea Mon Sep 17 00:00:00 2001 From: Yeicor Date: Sat, 8 Jan 2022 19:38:10 +0100 Subject: [PATCH] Heavily optimized STL importer --- examples/monkey_hat/main.go | 5 +-- go.mod | 1 + go.sum | 2 + obj/stl.go | 80 ++++++++++++++++++++++++++----------- sdf/voxel.go | 16 ++++---- 5 files changed, 70 insertions(+), 34 deletions(-) diff --git a/examples/monkey_hat/main.go b/examples/monkey_hat/main.go index 590b9f6d1..0b23d2be4 100644 --- a/examples/monkey_hat/main.go +++ b/examples/monkey_hat/main.go @@ -28,11 +28,10 @@ func monkeyWithHat() sdf.SDF3 { } } // - Create the SDF from the mesh (a modified Suzanne from Blender with 366 faces) - monkeyImported, err := obj.ImportSTL(file) + monkeyImported, err := obj.ImportSTL(file, 20, 3, 5) if err != nil { panic(err) } - //monkeyImported = sdf.NewVoxelSDF(monkeyImported, 64) // HAT hatHeight := 0.5 @@ -55,7 +54,7 @@ func monkeyWithHat() sdf.SDF3 { // It also smooths the mesh a little using trilinear interpolation. // It is actually slower for this mesh (unless meshCells <<< renderer's meshCells), but should be faster for // more complex meshes (with more triangles) or SDF3 hierarchies that take longer to evaluate. - monkeyHat = sdf.NewVoxelSDF(monkeyHat, 64, nil) // Use 32 for harder smoothing demo + monkeyHat = sdf.NewVoxelSDF3(monkeyHat, 64, nil) // Use 32 for harder smoothing demo return monkeyHat } diff --git a/go.mod b/go.mod index 4520fdd0c..95cebcb8e 100644 --- a/go.mod +++ b/go.mod @@ -4,6 +4,7 @@ go 1.13 require ( github.com/ajstarks/svgo v0.0.0-20200725142600-7a3c8b57fecb + github.com/dhconnelly/rtreego v1.1.0 github.com/golang/freetype v0.0.0-20170609003504-e2365dfdc4a0 github.com/hschendel/stl v1.0.4 github.com/llgcode/draw2d v0.0.0-20200930101115-bfaf5d914d1e diff --git a/go.sum b/go.sum index 2b9d1924e..41bf66bb8 100644 --- a/go.sum +++ b/go.sum @@ -7,6 +7,8 @@ github.com/ajstarks/svgo v0.0.0-20200725142600-7a3c8b57fecb/go.mod h1:K08gAheRH3 github.com/boombuler/barcode v1.0.0/go.mod h1:paBWMcWSl3LHKBqUq+rly7CNSldXjb2rDl3JlRe0mD8= github.com/davecgh/go-spew v1.1.0 h1:ZDRjVQ15GmhC3fiQ8ni8+OwkZQO4DARzQgrnXU1Liz8= github.com/davecgh/go-spew v1.1.0/go.mod h1:J7Y8YcW2NihsgmVo/mv3lAwl/skON4iLHjSsI+c5H38= +github.com/dhconnelly/rtreego v1.1.0 h1:ejMaqN03N1s6Bdg6peGkNgBnYYSBHzcK8yhSPCB+rHE= +github.com/dhconnelly/rtreego v1.1.0/go.mod h1:SDozu0Fjy17XH1svEXJgdYq8Tah6Zjfa/4Q33Z80+KM= github.com/fogleman/gg v1.2.1-0.20190220221249-0403632d5b90/go.mod h1:R/bRT+9gY/C5z7JzPU0zXsXHKM4/ayA+zqcVNZzPa1k= github.com/fogleman/gg v1.3.0/go.mod h1:R/bRT+9gY/C5z7JzPU0zXsXHKM4/ayA+zqcVNZzPa1k= github.com/go-fonts/dejavu v0.1.0/go.mod h1:4Wt4I4OU2Nq9asgDCteaAaWZOV24E+0/Pwo0gppep4g= diff --git a/obj/stl.go b/obj/stl.go index e808bb9db..3cfde91da 100644 --- a/obj/stl.go +++ b/obj/stl.go @@ -11,6 +11,7 @@ package obj import ( "github.com/deadsy/sdfx/render" "github.com/deadsy/sdfx/sdf" + "github.com/dhconnelly/rtreego" "github.com/hschendel/stl" "io" "math" @@ -19,22 +20,21 @@ import ( //----------------------------------------------------------------------------- type triMeshSdf struct { - tris []*render.Triangle3 - bb sdf.Box3 + rtree *rtreego.Rtree + numNeighbors int + bb sdf.Box3 } -const stlEpsilon = 1e-12 +const stlEpsilon = 1e-1 func (t *triMeshSdf) Evaluate(p sdf.V3) float64 { - if !t.bb.Contains(p) { // Fast exit - // Length to surface is at least distance to bounding box - return t.bb.Include(p).Size().Sub(t.bb.Size()).Length() - } // Check all triangle distances signedDistanceResult := 1. closestTriangle := math.MaxFloat64 - for _, triangle := range t.tris { - // TODO: Find a way to quickly skip this triangle (or a way to iterate a subset of triangles) + // Quickly skip checking most triangles by only checking the N closest neighbours (AABB based) + neighbors := t.rtree.NearestNeighbors(t.numNeighbors, stlToPoint(p)) + for _, neighbor := range neighbors { + triangle := neighbor.(*stlTriangle).Triangle3 testPointToTriangle := p.Sub(triangle.V[0]) triNormal := triangle.Normal() signedDistanceToTriPlane := triNormal.Dot(testPointToTriangle) @@ -52,15 +52,22 @@ func (t *triMeshSdf) BoundingBox() sdf.Box3 { return t.bb } -// ImportTriMesh converts a triangle-based mesh into a SDF3 surface. +// ImportTriMesh converts a triangle-based mesh into a SDF3 surface. minChildren and maxChildren are parameters that can +// affect the performance of the internal data structure (3 and 5 are a good default; maxChildren >= minChildren > 0). +// +// WARNING: Setting a low numNeighbors will consider many fewer triangles for each evaluated point, greatly speeding up +// the algorithm. However, if the count of triangles is too low artifacts will appear on the surface (triangle +// continuations). Setting this value to MaxInt is extremely slow but will provide correct results, so choose a value +// that works for your model. // -// It is recommended to cache its values at setup time by using sdf.VoxelSdf. +// It is recommended to cache (and/or smooth) its values by using sdf.VoxelSdf3. // // WARNING: It will only work on non-intersecting closed-surface(s) meshes. // NOTE: Fix using blender for intersecting surfaces: Edit mode > P > By loose parts > Add boolean modifier to join them -func ImportTriMesh(tris []*render.Triangle3) sdf.SDF3 { +func ImportTriMesh(tris chan *render.Triangle3, numNeighbors, minChildren, maxChildren int) sdf.SDF3 { m := &triMeshSdf{ - tris: tris, + rtree: nil, + numNeighbors: numNeighbors, bb: sdf.Box3{ Min: sdf.V3{X: math.MaxFloat64, Y: math.MaxFloat64, Z: math.MaxFloat64}, Max: sdf.V3{X: -math.MaxFloat64, Y: -math.MaxFloat64, Z: -math.MaxFloat64}, @@ -68,7 +75,9 @@ func ImportTriMesh(tris []*render.Triangle3) sdf.SDF3 { } // Compute the bounding box - for _, triangle := range tris { + bulkLoad := make([]rtreego.Spatial, 0) + for triangle := range tris { + bulkLoad = append(bulkLoad, &stlTriangle{Triangle3: triangle}) for _, vertex := range triangle.V { m.bb = m.bb.Include(vertex) } @@ -77,6 +86,7 @@ func ImportTriMesh(tris []*render.Triangle3) sdf.SDF3 { m.bb = sdf.Box3{} // Empty box centered at {0,0,0} } //m.bb = m.bb.ScaleAboutCenter(1 + 1e-12) // Avoids missing faces due to inaccurate math operations. + m.rtree = rtreego.NewTree(3, minChildren, maxChildren, bulkLoad...) return m } @@ -158,19 +168,43 @@ func stlPlaneProject(anyPoint, normal, testPoint sdf.V3) sdf.V3 { //----------------------------------------------------------------------------- +type stlTriangle struct { + *render.Triangle3 +} + +func (s *stlTriangle) Bounds() *rtreego.Rect { + bounds := sdf.Box3{Min: s.V[0], Max: s.V[0]} + bounds = bounds.Include(s.V[1]) + bounds = bounds.Include(s.V[2]) + points, err := rtreego.NewRectFromPoints(stlToPoint(bounds.Min), stlToPoint(bounds.Max)) + if err != nil { + panic(err) // Implementation error + } + return points +} + +func stlToPoint(v3 sdf.V3) rtreego.Point { + return rtreego.Point{v3.X, v3.Y, v3.Z} +} + +//----------------------------------------------------------------------------- + // ImportSTL converts an STL model into a SDF3 surface. See ImportTriMesh. -func ImportSTL(reader io.ReadSeeker) (sdf.SDF3, error) { +func ImportSTL(reader io.ReadSeeker, numNeighbors, minChildren, maxChildren int) (sdf.SDF3, error) { mesh, err := stl.ReadAll(reader) if err != nil { return nil, err } - tris := make([]*render.Triangle3, 0) // Buffer some triangles and send in batches if scheduler prefers it - for _, triangle := range mesh.Triangles { - tri := &render.Triangle3{} - for i, vertex := range triangle.Vertices { - tri.V[i] = sdf.V3{X: float64(vertex[0]), Y: float64(vertex[1]), Z: float64(vertex[2])} + tris := make(chan *render.Triangle3, 128) // Buffer some triangles and send in batches if scheduler prefers it + go func() { + for _, triangle := range mesh.Triangles { + tri := &render.Triangle3{} + for i, vertex := range triangle.Vertices { + tri.V[i] = sdf.V3{X: float64(vertex[0]), Y: float64(vertex[1]), Z: float64(vertex[2])} + } + tris <- tri } - tris = append(tris, tri) - } - return ImportTriMesh(tris), nil + close(tris) + }() + return ImportTriMesh(tris, numNeighbors, minChildren, maxChildren), nil } diff --git a/sdf/voxel.go b/sdf/voxel.go index 0c20659a8..2601ccd53 100644 --- a/sdf/voxel.go +++ b/sdf/voxel.go @@ -1,14 +1,14 @@ //----------------------------------------------------------------------------- /* -Voxel-based cache to remove deep SDF3 hierarchies at setup and speed up evaluation +Voxel-based cache/smoothing to remove deep SDF2/SDF3 hierarchies and speed up evaluation */ //----------------------------------------------------------------------------- package sdf -// VoxelSdf is the SDF that represents a pre-computed voxel-based SDF3. +// VoxelSdf3 is the SDF that represents a pre-computed voxel-based SDF3. //It can be used as a cache, or for smoothing. // // CACHE: @@ -18,7 +18,7 @@ package sdf // It performs trilinear mapping for inner values and may be used as a cache for any other SDF, losing some accuracy. // // WARNING: It may lose sharp features, even if meshCells is high. -type VoxelSdf struct { +type VoxelSdf3 struct { // voxelCorners are the values of this SDF in each voxel corner voxelCorners map[V3i]float64 // TODO: Octree + k-d tree to simplify/reduce memory consumption + speed-up access? // bb is the bounding box. @@ -27,8 +27,8 @@ type VoxelSdf struct { numVoxels V3i } -// NewVoxelSDF see VoxelSdf. This populates the whole cache from the given SDF. The progress listener may be nil. -func NewVoxelSDF(s SDF3, meshCells int, progress chan float64) SDF3 { +// NewVoxelSDF3 see VoxelSdf3. This populates the whole cache from the given SDF. The progress listener may be nil. +func NewVoxelSDF3(s SDF3, meshCells int, progress chan float64) SDF3 { bb := s.BoundingBox() // TODO: Use default code to avoid duplication bbSize := bb.Size() resolution := bbSize.MaxComponent() / float64(meshCells) @@ -48,14 +48,14 @@ func NewVoxelSDF(s SDF3, meshCells int, progress chan float64) SDF3 { } } - return &VoxelSdf{ + return &VoxelSdf3{ voxelCorners: voxelCorners, bb: bb, numVoxels: cells, } } -func (m *VoxelSdf) Evaluate(p V3) float64 { +func (m *VoxelSdf3) Evaluate(p V3) float64 { // Find the voxel's {0,0,0} corner quickly and compute p's displacement voxelSize := m.bb.Size().Div(m.numVoxels.ToV3()) voxelStartIndex := p.Sub(m.bb.Min).Div(voxelSize).ToV3i() @@ -84,6 +84,6 @@ func (m *VoxelSdf) Evaluate(p V3) float64 { return c } -func (m *VoxelSdf) BoundingBox() Box3 { +func (m *VoxelSdf3) BoundingBox() Box3 { return m.bb }