merge heads

5 files changed, 89 insertions, 135 deletions

diff --git a/camera.py b/camera.py
index 59c47eb..005e2cd 100755
--- a/camera.py
+++ b/camera.py
@@ -72,11 +72,17 @@ def get_rays(position, size = (800, 600), film_size = (0.035, 0.024), focal_leng
     return grid, focal_point-grid
 
 class Camera(Thread):
-    def __init__(self, geometry, module, context, lock=None, size=(800,600)):
+    def __init__(self, geometry, module, context, lock=None, size=(800,600),
+                 spnav=False):
         Thread.__init__(self)
         self.geometry = geometry
         self.module = module
         self.context = context
+        self.spnav = spnav
+
+        if spnav:
+            import spnav as spnav_module
+            self.spnav_module = spnav_module
 
         if lock is None:
             self.lock = Lock()
@@ -130,6 +136,7 @@ class Camera(Thread):
     def initialize_render(self):
         with self.lock:
             self.context.push()
+            print 'beginning initialize_render()'
             self.rng_states_gpu = cuda.mem_alloc(self.width*self.height*sizeof('curandStateXORWOW', '#include <curand_kernel.h>'))
             self.init_rng_kernel(np.int32(self.width*self.height), self.rng_states_gpu, np.int32(0), np.int32(0), block=(self.nblocks,1,1), grid=(self.width*self.height//self.nblocks+1,1))
             self.xyz_lookup1_gpu = gpuarray.zeros(len(self.geometry.mesh.triangles), dtype=gpuarray.vec.float3)
@@ -222,7 +229,7 @@ class Camera(Thread):
 
         self.update()
 
-    def rotate_around_point(self, phi, n, point):
+    def rotate_around_point(self, phi, n, point, redraw=True):
         with self.lock:
             self.context.push()
             self.rotate_around_point_kernel(np.int32(self.origins_gpu.size), self.origins_gpu, np.float32(phi), gpuarray.vec.make_float3(*n), gpuarray.vec.make_float3(*point), block=(self.nblocks,1,1), grid=(self.origins_gpu.size//self.nblocks+1,1))
@@ -232,12 +239,13 @@ class Camera(Thread):
         self.axis1 = rotate(self.axis1, phi, n)
         self.axis2 = rotate(self.axis2, phi, n)
 
-        if self.render:
-            self.clear_image()
+        if redraw:
+            if self.render:
+                self.clear_image()
 
-        self.update()
+            self.update()
 
-    def translate(self, v):
+    def translate(self, v, redraw=True):
         with self.lock:
             self.context.push()
             self.translate_kernel(np.int32(self.pixels_gpu.size), self.origins_gpu, gpuarray.vec.make_float3(*v), block=(self.nblocks,1,1), grid=(self.pixels_gpu.size//self.nblocks,1))
@@ -245,14 +253,15 @@ class Camera(Thread):
             self.point += v
             self.context.pop()
 
-        if self.render:
-            self.clear_image()
+        if redraw:
+            if self.render:
+                self.clear_image()
 
-        self.update()
+            self.update()
 
     def update(self):
         if self.render:
-            while self.nlookup_calls < 100:
+            while self.nlookup_calls < 10:
                 self.update_xyz_lookup(self.source_position)
             self.update_image()
             self.process_image()
@@ -273,6 +282,9 @@ class Camera(Thread):
             self.movie_index += 1
 
     def run(self):
+        if self.spnav:
+            self.spnav_module.spnav_open()
+
         self.update()
         
         done = False
@@ -287,6 +299,58 @@ class Camera(Thread):
 
             if self.render and not clicked and not pygame.event.peek(KEYDOWN):
                 self.update()
+            
+            # Space Navigator controls
+            if self.spnav:
+                spnav_event = self.spnav_module.spnav_poll_event()
+                if spnav_event:
+                    if spnav_event.type == self.spnav_module.SPNAV_EVENT_MOTION:
+                        if shift:
+                            accelerate_factor = 2.0
+                        else:
+                            accelerate_factor = 1.0
+                        #print 'raw:', spnav_event
+                        v1 = self.axis1
+                        v2 = self.axis2
+                        v3 = np.cross(self.axis1,self.axis2)
+
+                        v = -v1*spnav_event.motion.x + v2*spnav_event.motion.y \
+                            + v3*spnav_event.motion.z
+                        v *= self.scale / 10000.0 * accelerate_factor
+
+                        #print 'translate:', v
+                        self.translate(v, redraw=False)                    
+
+                        axis = v1*spnav_event.motion.rx + -v2*spnav_event.motion.ry \
+                               + -v3*spnav_event.motion.rz
+
+                        # Zero all but non-max values
+                        #axis[~(np.abs(axis) == np.abs(axis).max())] = 0
+                        #axis[np.abs(axis) < 15] = 0
+                        if (axis != 0).any():
+                            axis = axis.astype(float)
+                            length = np.linalg.norm(axis)
+                            angle = length * 0.0001 * accelerate_factor
+                            axis /= length
+                            #print 'rotate:', angle, axis
+                            self.rotate_around_point(angle, axis, self.point, redraw=False)
+
+                        if self.render:
+                            self.clear_image()
+
+                        self.update()
+                        self.spnav_module.spnav_remove_events(self.spnav_module.SPNAV_EVENT_MOTION)
+
+                    elif spnav_event.type == self.spnav_module.SPNAV_EVENT_BUTTON:
+                        if spnav_event.button.bnum == 0 and spnav_event.button.press:
+                            if not hasattr(self, 'rng_states_gpu'):
+                                self.initialize_render()
+
+                            self.render = not self.render
+                            self.clear_image()
+                            self.update()
+                            self.spnav_module.spnav_remove_events(self.spnav_module.SPNAV_EVENT_ANY)
+                        
 
             for event in pygame.event.get():
                 if event.type == MOUSEBUTTONDOWN:
@@ -436,6 +500,9 @@ class Camera(Thread):
             encode_movie(self.movie_dir)
 
         pygame.display.quit()
+        if self.spnav:
+            self.spnav_module.spnav_close()
+
 
 if __name__ == '__main__':
     import optparse
@@ -452,11 +519,13 @@ if __name__ == '__main__':
 
     parser = optparse.OptionParser('%prog filename.stl')
     parser.add_option('-b', '--bits', type='int', dest='bits',
-                      help='bits for z-ordering space axes', default=8)
+                      help='bits for z-ordering space axes', default=10)
     #parser.add_option('-l', action='store_true', dest='load_bvh',
     #                  help='load bounding volumes', default=False)
     parser.add_option('-r', '--resolution', dest='resolution',
                       help='specify resolution', default='800,600')
+    parser.add_option('--spnav', action='store_true', dest='spnav',
+                      help='activate Space Navigator support', default=False)
     options, args = parser.parse_args()
 
     if len(args) < 1:
@@ -493,5 +562,6 @@ if __name__ == '__main__':
     gpu = GPU()
     gpu.load_geometry(geometry)
     gpu.context.pop()
-    camera = Camera(geometry, gpu.module, gpu.context, lock, size=size)
+    camera = Camera(geometry, gpu.module, gpu.context, lock, size=size,
+                    spnav=options.spnav)
     camera.start()
diff --git a/geometry.py b/geometry.py
index b513112..1b85c51 100644
--- a/geometry.py
+++ b/geometry.py
@@ -339,120 +339,3 @@ class Geometry(object):
             if unique_bit_shifted_zvalues.size == 1:
                 break
 
-    def load(self, module):
-        """
-        Load the bounding volume hierarchy onto the GPU module,
-        bind it to the appropriate textures, and return a list
-        of the texture references.
-        """
-
-        if not hasattr(self, 'mesh'):
-            raise RuntimeError('cannot load before a call to build().')
-
-        set_wavelength_range = module.get_function('set_wavelength_range')
-
-        set_wavelength_range(np.float32(standard_wavelengths[0]), np.float32(standard_wavelengths[-1]), np.float32(standard_wavelengths[1]-standard_wavelengths[0]), np.uint32(standard_wavelengths.size), block=(1,1,1), grid=(1,1))
-
-        set_material = module.get_function('set_material')
-        for i, material in enumerate(self.unique_materials):
-            if material is None:
-                if color is False:
-                    raise Exception('one or more triangles is missing a material.')
-                continue
-
-            refractive_index = np.interp(standard_wavelengths, material.refractive_index[:,0], material.refractive_index[:,1]).astype(np.float32)
-            absorption_length = np.interp(standard_wavelengths, material.absorption_length[:,0], material.absorption_length[:,1]).astype(np.float32)
-            scattering_length = np.interp(standard_wavelengths, material.scattering_length[:,0], material.scattering_length[:,1]).astype(np.float32) 
-
-            material.refractive_index_gpu = cuda.to_device(refractive_index)
-            material.absorption_length_gpu = cuda.to_device(absorption_length)
-            material.scattering_length_gpu = cuda.to_device(scattering_length)
-
-            set_material(np.int32(i), material.refractive_index_gpu, material.absorption_length_gpu, material.scattering_length_gpu, block=(1,1,1), grid=(1,1))
-
-        set_surface = module.get_function('set_surface')
-        for i, surface in enumerate(self.unique_surfaces):
-            if surface is None:
-                continue
-
-            detect = np.interp(standard_wavelengths, surface.detect[:,0], surface.detect[:,1]).astype(np.float32)
-            absorb = np.interp(standard_wavelengths, surface.absorb[:,0], surface.absorb[:,1]).astype(np.float32)
-            reflect_diffuse = np.interp(standard_wavelengths, surface.reflect_diffuse[:,0], surface.reflect_diffuse[:,1]).astype(np.float32)
-            reflect_specular = np.interp(standard_wavelengths, surface.reflect_specular[:,0], surface.reflect_specular[:,1]).astype(np.float32)
-
-            surface.detect_gpu = cuda.to_device(detect)
-            surface.absorb_gpu = cuda.to_device(absorb)
-            surface.reflect_diffuse_gpu = cuda.to_device(reflect_diffuse)
-            surface.reflect_specular_gpu = cuda.to_device(reflect_specular)
-
-            set_surface(np.int32(i), surface.detect_gpu, surface.absorb_gpu, surface.reflect_diffuse_gpu, surface.reflect_specular_gpu, block=(1,1,1), grid=(1,1))
-
-        vertices = np.empty(len(self.mesh.vertices), dtype=gpuarray.vec.float3)
-        vertices['x'] = self.mesh.vertices[:,0]
-        vertices['y'] = self.mesh.vertices[:,1]
-        vertices['z'] = self.mesh.vertices[:,2]
-
-        triangles = \
-            np.empty(len(self.mesh.triangles), dtype=gpuarray.vec.uint4)
-        triangles['x'] = self.mesh.triangles[:,0]
-        triangles['y'] = self.mesh.triangles[:,1]
-        triangles['z'] = self.mesh.triangles[:,2]
-        triangles['w'] = ((self.material1_index & 0xff) << 24) | ((self.material2_index & 0xff) << 16) | ((self.surface_index & 0xff) << 8)
-
-        lower_bounds = np.empty(self.lower_bounds.shape[0], dtype=gpuarray.vec.float4)
-        lower_bounds['x'] = self.lower_bounds[:,0]
-        lower_bounds['y'] = self.lower_bounds[:,1]
-        lower_bounds['z'] = self.lower_bounds[:,2]
-
-        upper_bounds = np.empty(self.upper_bounds.shape[0], dtype=gpuarray.vec.float4)
-        upper_bounds['x'] = self.upper_bounds[:,0]
-        upper_bounds['y'] = self.upper_bounds[:,1]
-        upper_bounds['z'] = self.upper_bounds[:,2]
-
-        self.vertices_gpu = cuda.to_device(vertices)
-        self.triangles_gpu = cuda.to_device(triangles)
-        self.colors_gpu = gpuarray.to_gpu(self.colors.astype(np.uint32))
-        self.lower_bounds_gpu = cuda.to_device(lower_bounds)
-        self.upper_bounds_gpu = cuda.to_device(upper_bounds)
-        self.node_map_gpu = cuda.to_device(self.node_map)
-        self.node_length_gpu = cuda.to_device(self.node_length)
-
-        def format_size(size):
-            if size < 1e3:
-                return '%.1f%s' % (size, ' ')
-            elif size < 1e6:
-                return '%.1f%s' % (size/1e3, 'K')
-            elif size < 1e9:
-                return '%.1f%s' % (size/1e6, 'M')
-            else:
-                return '%.1f%s' % (size/1e9, 'G')
-
-        def format_array(name, array):
-            return '%-15s %6s %6s' % (name, format_size(len(array)), format_size(array.nbytes))
-
-        print 'Device Usage:'
-        print format_array('vertices', vertices)
-        print format_array('triangles', triangles)
-        print format_array('lower_bounds', self.lower_bounds)
-        print format_array('upper_bounds', self.upper_bounds)
-        print format_array('node_map', self.node_map)
-        print format_array('node_length', self.node_length)
-        print '%-15s %6s %6s' % ('total', '', format_size(vertices.nbytes + triangles.nbytes + self.lower_bounds.nbytes + self.upper_bounds.nbytes + self.node_map.nbytes + self.node_length.nbytes))
-
-        set_global_mesh_variables = module.get_function('set_global_mesh_variables')
-        set_global_mesh_variables(self.triangles_gpu, self.vertices_gpu, self.colors_gpu, np.uint32(self.node_map.size-1), np.uint32(self.first_node), block=(1,1,1), grid=(1,1))
-
-        self.lower_bounds_tex = module.get_texref('lower_bounds')
-        self.upper_bounds_tex = module.get_texref('upper_bounds')
-        self.node_map_tex = module.get_texref('node_map')
-        self.node_length_tex = module.get_texref('node_length')
-
-        self.lower_bounds_tex.set_address(self.lower_bounds_gpu, lower_bounds.nbytes)
-        self.upper_bounds_tex.set_address(self.upper_bounds_gpu, upper_bounds.nbytes)
-        self.node_map_tex.set_address(self.node_map_gpu, self.node_map.nbytes)
-        self.node_length_tex.set_address(self.node_length_gpu, self.node_length.nbytes)
-
-        self.lower_bounds_tex.set_format(cuda.array_format.FLOAT, 4)
-        self.upper_bounds_tex.set_format(cuda.array_format.FLOAT, 4)
-        self.node_map_tex.set_format(cuda.array_format.UNSIGNED_INT32, 1)
-        self.node_length_tex.set_format(cuda.array_format.UNSIGNED_INT32, 1)
diff --git a/gpu.py b/gpu.py
index 11a83ee..e5040a8 100644
--- a/gpu.py
+++ b/gpu.py
@@ -74,6 +74,7 @@ class GPU(object):
 
     def load_geometry(self, geometry):
         if not hasattr(geometry, 'mesh'):
+            print 'WARNING: building geometry with 8-bits'
             geometry.build(bits=8)
 
         self.geo_funcs.set_wavelength_range(np.float32(standard_wavelengths[0]), np.float32(standard_wavelengths[-1]), np.float32(standard_wavelengths[1]-standard_wavelengths[0]), np.uint32(standard_wavelengths.size), block=(1,1,1), grid=(1,1))
diff --git a/optics.py b/optics.py
index dea5a20..3580dc4 100644
--- a/optics.py
+++ b/optics.py
@@ -3,8 +3,8 @@ from geometry import Material, Surface
 
 vacuum = Material('vacuum')
 vacuum.set('refractive_index', 1.0)
-vacuum.set('absorption_length', np.inf)
-vacuum.set('scattering_length', np.inf)
+vacuum.set('absorption_length', 1e6)
+vacuum.set('scattering_length', 1e6)
 
 lambertian_surface = Surface('lambertian_surface')
 lambertian_surface.set('reflect_diffuse', 1)
@@ -153,7 +153,7 @@ glass = Material('glass')
 glass.set('refractive_index', 1.49)
 glass.absorption_length = \
     np.array([(200, 0.1e-6), (300, 0.1e-6), (330, 1.0), (500, 2.0), (600, 1.0), (770, 0.5), (800, 0.1e-6)])
-glass.set('scattering_length', np.inf)
+glass.set('scattering_length', 1e6)
 
 # From SNO+ database
 acrylic_sno = Material('acrylic_sno')
diff --git a/sim.py b/sim.py
index de219c8..6432e7d 100755
--- a/sim.py
+++ b/sim.py
@@ -39,7 +39,7 @@ if __name__ == '__main__':
     parser.add_option('--particle', type='string', dest='particle', default='e-')
     parser.add_option('--energy', type='float', dest='energy', default=100.0)
     parser.add_option('--pos', type='string', dest='pos', default='(0,0,0)')
-    parser.add_option('--dir', type='string', dest='pos', default='(1,0,0)')
+    parser.add_option('--dir', type='string', dest='dir', default='(1,0,0)')
 
     options, args = parser.parse_args()
     if len(args) != 1:
@@ -53,7 +53,7 @@ if __name__ == '__main__':
         sys.exit(1)
 
     position = np.array(eval(options.pos), dtype=float)
-    direction = np.array(eval(options.pos), dtype=float)
+    direction = np.array(eval(options.dir), dtype=float)
     detector = detectors.find(options.detector)
 
     print >>sys.stderr, 'Creating detector...'
@@ -102,7 +102,7 @@ if __name__ == '__main__':
                                              position=position,
                                              direction=direction)
         nphotons += len(photons['pos'])
-
+ 
         # this will stop and wait for event to finish
         hit_times = gpu_worker.get_hits()