1 files changed, 75 insertions, 37 deletions
diff --git a/geometry.py b/geometry.py
index d365141..a65dbc4 100644
--- a/geometry.py
+++ b/geometry.py
@@ -1,63 +1,101 @@
 import numpy as np
-from bvh import *
+from itertools import chain
 
-class Geometry(object):
-    """
-    Geometry object.
-    """
+def compress(data, selectors):
+    return (d for d, s in zip(data, selectors) if s)
 
-    def __init__(self):
-        self.materials = []
-        self.mesh = np.empty(0)
-        self.inside = np.empty(0)
-        self.outside = np.empty(0)
+class Leaf(object):
+    def __init__(self, mesh, mesh_idx, zvalue=None):
+        self.zvalue = zvalue
+        self.mesh_idx = mesh_idx
+
+        self.lower_bound = np.array([np.min(mesh[:,:,0]), np.min(mesh[:,:,1]), np.min(mesh[:,:,2])])
+        self.upper_bound = np.array([np.max(mesh[:,:,0]), np.max(mesh[:,:,1]), np.max(mesh[:,:,2])])
+
+        self.size = mesh.shape[0]
+
+class Node(object):
+    def __init__(self, children, zvalue=None):
+        self.zvalue = zvalue
+
+        lower_pts = np.array([child.lower_bound for child in children])
+        upper_pts = np.array([child.upper_bound for child in children])
+
+        self.lower_bound = np.array([np.min(lower_pts[:,0]), np.min(lower_pts[:,1]), np.min(lower_pts[:,2])])
+        self.upper_bound = np.array([np.max(upper_pts[:,0]), np.max(upper_pts[:,1]), np.max(upper_pts[:,2])])
+
+        self.size = len(children)
+
+        self.children = children
 
-    def add_mesh(self, mesh, inside, outside):
-        """
-        Add a triangle mesh.
+def interleave(*args):
+    return int("".join(chain(*zip(*[bin(x)[2:].zfill(x.nbytes*8) for x in args]))), 2)
 
-        Args:
-            - mesh: array of shape (n,3,3)
-                Array of vertices for a closed mesh with n triangles.
+def morton_order(mesh, bits=3):
+    lower_bound = np.array([np.min(mesh[:,:,0]), np.min(mesh[:,:,1]), np.min(mesh[:,:,2])])
+    upper_bound = np.array([np.max(mesh[:,:,0]), np.max(mesh[:,:,1]), np.max(mesh[:,:,2])])
 
-            - inside: Material
-                Material inside the mesh.
+    max_value = 2**bits - 1
+    
+    def quantize(x):
+        return np.uint64((x-lower_bound)*max_value/(upper_bound-lower_bound))
 
-            - outside: Material
-                Material outside the mesh.
-         """
+    zvalues = np.empty(mesh.shape[0], dtype=np.uint64)
+    for i, triangle in enumerate(mesh):
+        zvalues[i] = interleave(*quantize(np.mean(triangle, axis=0)))
 
+    return zvalues
+
+class Solid(object):
+    def __init__(self, mesh, inside, outside):
         if len(mesh.shape) != 3 or mesh.shape[1] != 3 or mesh.shape[2] != 3:
             raise Exception('shape mismatch')
 
-        if inside not in self.materials:
-            self.materials.append(inside)
-        if outside not in self.materials:
-            self.materials.append(outside)
+        self.mesh = mesh
+        self.inside = inside
+        self.outside = outside
+
+    def __len__(self):
+        return self.mesh.shape[0]
+
+class Geometry(object):
+    """
+    Geometry object.
+    """
+
+    def __init__(self):
+        self.solids = []
+        self.materials = []
+
+    def add_solid(self, solid):
+        self.solids.append(solid)
 
-        self.inside.resize(self.inside.size+mesh.shape[0])
-        self.inside[-mesh.shape[0]:] = np.tile(self.materials.index(inside), mesh.shape[0])
+        if solid.inside not in self.materials:
+            self.materials.append(solid.inside)
 
-        self.outside.resize(self.outside.size+mesh.shape[0])
-        self.outside[-mesh.shape[0]:] = np.tile(self.materials.index(outside), mesh.shape[0])
+        if solid.outside not in self.materials:
+            self.materials.append(solid.outside)
 
-        self.mesh.resize((self.mesh.shape[0]+mesh.shape[0],3,3))
-        self.mesh[-mesh.shape[0]:] = mesh
+    def build(self, bits=3):
+        self.mesh = np.concatenate([solid.mesh for solid in self.solids])
+        self.solid_index = np.concatenate([np.tile(self.solids.index(solid), len(solid)) for solid in self.solids])
+        self.inside_index = np.concatenate([np.tile(self.materials.index(solid.outside), len(solid)) for solid in self.solids])
+        self.outside_index = np.concatenate([np.tile(self.materials.index(solid.outside), len(solid)) for solid in self.solids])
 
-    def build(self, bits=16):
-        order, zvalues = morton_order(self.mesh, bits)
+        zvalues = morton_order(self.mesh, bits)
 
+        order = np.array(zip(*sorted(zip(zvalues, range(zvalues.size))))[-1])
 
         zvalues = zvalues[order]
         self.mesh = self.mesh[order]
-        self.inside = self.inside[order]
-        self.outside = self.outside[order]
+        self.solid_index = self.solid_index[order]
+        self.inside_index = self.inside_index[order]
+        self.outside_index = self.outside_index[order]
 
         leafs = []
-
         for z in sorted(set(zvalues)):
             mask = (zvalues == z)
-            leafs.append(Leaf(self.mesh[mask], 3*np.where(mask)[0][0], z))
+            leafs.append(Leaf(self.mesh[mask], np.where(mask)[0][0], z))
 
         layers = []
         layers.append(leafs)