From 1a6dc30108d3e78f72f773ec025fc98e246f68f5 Mon Sep 17 00:00:00 2001 From: Anthony LaTorre Date: Sat, 30 Jul 2011 16:10:27 -0400 Subject: when throwing photons from the light source out onto the scene, photons are now thrown randomly across each triangle instead of only at the center of each triangle. all of the rendering kernels have been rewritten so that they operate additively; for example, you may now throw photons from the light source onto the scene, render from the camera to the scene, then throw more photons and render again. --- camera.py | 175 ++++++++--------- color/Illuminanta.csv | 531 ++++++++++++++++++++++++++++++++++++++++++++++++++ color/chromaticity.py | 18 +- color/scvle_1.csv | 401 ++++++++++++++++++++++++++++++++++++++ color/vl1924e_1.csv | 471 ++++++++++++++++++++++++++++++++++++++++++++ detectors/__init__.py | 4 + geometry.py | 14 +- solids/__init__.py | 6 + src/kernel.cu | 238 ++++++++++------------ src/photon.h | 8 +- 10 files changed, 1607 insertions(+), 259 deletions(-) create mode 100644 color/Illuminanta.csv create mode 100644 color/scvle_1.csv create mode 100644 color/vl1924e_1.csv diff --git a/camera.py b/camera.py index a786971..148ddf9 100644 --- a/camera.py +++ b/camera.py @@ -2,8 +2,10 @@ import numpy as np from itertools import product, count from threading import Thread, Lock import time +import datetime import os import sys +from color import map_wavelength from transform import rotate @@ -14,6 +16,14 @@ from pycuda import gpuarray from pycuda.characterize import sizeof import pycuda.driver as cuda +def timeit(func): + def f(*args, **kwargs): + t0 = time.time() + func(*args, **kwargs) + elapsed = time.time() - t0 + print '%s elapsed in %s().' % (datetime.timedelta(seconds=elapsed), func.__name__) + return f + def get_rays(position, size = (800, 600), film_size = (0.035, 0.024), focal_length=0.05): """ Generate ray positions and directions from a pinhole camera facing the negative y direction. @@ -63,8 +73,9 @@ class Camera(Thread): self.ray_trace_kernel = self.module.get_function('ray_trace') self.rotate_kernel = self.module.get_function('rotate') self.translate_kernel = self.module.get_function('translate') - self.build_rgb_lookup_kernel = self.module.get_function('build_rgb_lookup') - self.render_kernel = self.module.get_function('render') + self.update_xyz_lookup_kernel = self.module.get_function('update_xyz_lookup') + self.update_xyz_image_kernel = self.module.get_function('update_xyz_image') + self.process_image_kernel = self.module.get_function('process_image') self.init_rng_kernel = self.module.get_function('init_rng') self.context.pop() @@ -77,7 +88,7 @@ class Camera(Thread): self.axis1 = np.array([1,0,0], float) self.axis2 = np.array([0,0,1], float) - origins, directions = get_rays(self.point) + origins, directions = get_rays(self.point, self.size) with self.lock: self.context.push() @@ -88,77 +99,73 @@ class Camera(Thread): self.render = False - def build_rgb_lookup(self, source_position=(0,0,0)): + @timeit + def initialize_render(self): with self.lock: self.context.push() - print '\ninitializing random number states.' self.rng_states_gpu = cuda.mem_alloc(self.width*self.height*sizeof('curandStateXORWOW', '#include ')) 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) + self.xyz_lookup2_gpu = gpuarray.zeros(len(self.geometry.mesh.triangles), dtype=gpuarray.vec.float3) + self.image_gpu = gpuarray.zeros(self.width*self.height, dtype=gpuarray.vec.float3) self.context.synchronize() - print 'done.' - - self.source_positions_gpu = gpuarray.empty(len(geometry.mesh.triangles), dtype=gpuarray.vec.float3) - self.source_positions_gpu.fill(gpuarray.vec.make_float3(*source_position)) + self.context.pop() - source_directions = np.mean(self.geometry.mesh[:], axis=1) - source_position - self.source_directions_gpu = gpuarray.to_gpu(source_directions.astype(np.float32).view(gpuarray.vec.float3)) + self.source_position = self.point - self.rgb_lookup1_gpu = gpuarray.zeros(self.source_positions_gpu.size, dtype=gpuarray.vec.float3) - self.rgb_lookup2_gpu = gpuarray.zeros(self.source_positions_gpu.size, dtype=gpuarray.vec.float3) + self.nimages = 0 + self.nlookup_calls = 0 + self.max_steps = 10 - self.max_steps = 10 - rgb_runs = 100 - - print 'building rgb lookup.' - self.build_rgb_lookup_kernel(np.int32(self.source_positions_gpu.size), self.source_positions_gpu, self.source_directions_gpu, self.rgb_lookup1_gpu, self.rgb_lookup2_gpu, np.uint64(np.random.randint(np.iinfo(np.int64).max)), np.int32(rgb_runs), np.int32(self.max_steps), block=(self.nblocks,1,1), grid=(self.source_positions_gpu.size//self.nblocks+1,1)) - self.context.synchronize() - print 'done.' + def clear_xyz_lookup(self): + with self.lock: + self.context.push() + self.xyz_lookup1_gpu.fill(gpuarray.vec.make_float3(0.0,0.0,0.0)) + self.xyz_lookup2_gpu.fill(gpuarray.vec.make_float3(0.0,0.0,0.0)) + self.context.pop() - rgb_lookup1 = self.rgb_lookup1_gpu.get().view(np.float32) - rgb_lookup1 /= rgb_runs - rgb_lookup1[rgb_lookup1 > 1.0] = 1.0 - self.rgb_lookup1_gpu.set(rgb_lookup1.view(gpuarray.vec.float3)) + self.nlookup_calls = 0 - rgb_lookup2 = self.rgb_lookup2_gpu.get().view(np.float32) - rgb_lookup2 /= rgb_runs - rgb_lookup2[rgb_lookup2 > 1.0] = 1.0 - self.rgb_lookup2_gpu.set(rgb_lookup2.view(gpuarray.vec.float3)) + def update_xyz_lookup(self, source_position): + with self.lock: + self.context.push() + for i in range(self.xyz_lookup1_gpu.size//(self.width*self.height)+1): + self.update_xyz_lookup_kernel(np.int32(self.width*self.height), np.int32(self.xyz_lookup1_gpu.size), np.int32(i*self.width*self.height), gpuarray.vec.make_float3(*source_position), self.rng_states_gpu, np.float32(685.0), gpuarray.vec.make_float3(1.0,0.0,0.0), self.xyz_lookup1_gpu, self.xyz_lookup2_gpu, np.int32(self.max_steps), block=(self.nblocks,1,1), grid=(self.width*self.height//self.nblocks+1,1)) - self.nimages = 0 - self.rgb_image = np.zeros((self.pixels_gpu.size, 3), float) - self.r_gpu = gpuarray.zeros(self.pixels_gpu.size, dtype=np.float32) - self.g_gpu = gpuarray.zeros(self.pixels_gpu.size, dtype=np.float32) - self.b_gpu = gpuarray.zeros(self.pixels_gpu.size, dtype=np.float32) + for i in range(self.xyz_lookup1_gpu.size//(self.width*self.height)+1): + self.update_xyz_lookup_kernel(np.int32(self.width*self.height), np.int32(self.xyz_lookup1_gpu.size), np.int32(i*self.width*self.height), gpuarray.vec.make_float3(*source_position), self.rng_states_gpu, np.float32(545.0), gpuarray.vec.make_float3(0.0,1.0,0.0), self.xyz_lookup1_gpu, self.xyz_lookup2_gpu, np.int32(self.max_steps), block=(self.nblocks,1,1), grid=(self.width*self.height//self.nblocks+1,1)) + for i in range(self.xyz_lookup1_gpu.size//(self.width*self.height)+1): + self.update_xyz_lookup_kernel(np.int32(self.width*self.height), np.int32(self.xyz_lookup1_gpu.size), np.int32(i*self.width*self.height), gpuarray.vec.make_float3(*source_position), self.rng_states_gpu, np.float32(445.0), gpuarray.vec.make_float3(0.0,0.0,1.0), self.xyz_lookup1_gpu, self.xyz_lookup2_gpu, np.int32(self.max_steps), block=(self.nblocks,1,1), grid=(self.width*self.height//self.nblocks+1,1)) self.context.pop() + self.nlookup_calls += 1 + def clear_image(self): - try: - self.rgb_image.fill(0.0) - self.nimages = 0 - except AttributeError: - pass + with self.lock: + self.context.push() + self.image_gpu.fill(gpuarray.vec.make_float3(0.0,0.0,0.0)) + self.context.pop() - def acquire_image(self): - if not hasattr(self, 'rng_states_gpu'): - self.build_rgb_lookup() + self.nimages = 0 + def update_image(self): with self.lock: self.context.push() - self.render_kernel(np.int32(self.pixels_gpu.size), self.rng_states_gpu, self.origins_gpu, self.directions_gpu, self.r_gpu, self.g_gpu, self.b_gpu, self.rgb_lookup1_gpu, self.rgb_lookup2_gpu, np.int32(self.max_steps), block=(self.nblocks,1,1), grid=(self.pixels_gpu.size//self.nblocks+1,1)) - self.rgb_image[:,0] += self.r_gpu.get() - self.rgb_image[:,1] += self.g_gpu.get() - self.rgb_image[:,2] += self.b_gpu.get() - self.nimages += 1 - self.context.pop() + self.update_xyz_image_kernel(np.int32(self.width*self.height), self.rng_states_gpu, self.origins_gpu, self.directions_gpu, np.float32(685.0), gpuarray.vec.make_float3(1.0,0.0,0.0), self.xyz_lookup1_gpu, self.xyz_lookup2_gpu, self.image_gpu, np.int32(self.nlookup_calls), np.int32(self.max_steps), block=(self.nblocks,1,1), grid=(self.width*self.height//self.nblocks+1,1)) - def render_image(self): - scaled_rgb_image = ((self.rgb_image/self.nimages)*255).astype(np.int32) + self.update_xyz_image_kernel(np.int32(self.width*self.height), self.rng_states_gpu, self.origins_gpu, self.directions_gpu, np.float32(545.0), gpuarray.vec.make_float3(0.0,1.0,0.0), self.xyz_lookup1_gpu, self.xyz_lookup2_gpu, self.image_gpu, np.int32(self.nlookup_calls), np.int32(self.max_steps), block=(self.nblocks,1,1), grid=(self.width*self.height//self.nblocks+1,1)) - image = scaled_rgb_image[:,0] << 16 | scaled_rgb_image[:,1] << 8 | scaled_rgb_image[:,2] + self.update_xyz_image_kernel(np.int32(self.width*self.height), self.rng_states_gpu, self.origins_gpu, self.directions_gpu, np.float32(445.0), gpuarray.vec.make_float3(0.0,0.0,1.0), self.xyz_lookup1_gpu, self.xyz_lookup2_gpu, self.image_gpu, np.int32(self.nlookup_calls), np.int32(self.max_steps), block=(self.nblocks,1,1), grid=(self.width*self.height//self.nblocks+1,1)) + self.context.pop() - pygame.surfarray.blit_array(self.screen, image.reshape(self.size)) - pygame.display.flip() + self.nimages += 1 + + def process_image(self): + with self.lock: + self.context.push() + self.process_image_kernel(np.int32(self.width*self.height), self.image_gpu, self.pixels_gpu, np.int32(self.nimages), block=(self.nblocks,1,1), grid=((self.width*self.height)//self.nblocks+1,1)) + self.context.pop() def screenshot(self, dir=''): root, ext = 'screenshot', 'png' @@ -183,6 +190,11 @@ class Camera(Thread): self.axis2 = rotate(self.axis2, phi, n) self.context.pop() + if self.render: + self.clear_image() + + self.update() + def translate(self, v): with self.lock: self.context.push() @@ -191,22 +203,25 @@ class Camera(Thread): self.point += v self.context.pop() + if self.render: + self.clear_image() + + self.update() + def update(self): if self.render: - self.acquire_image() - self.render_image() - return + while self.nlookup_calls < 10: + self.update_xyz_lookup(self.source_position) + self.update_image() + self.process_image() + else: + with self.lock: + self.context.push() + self.ray_trace_kernel(np.int32(self.pixels_gpu.size), self.origins_gpu, self.directions_gpu, self.pixels_gpu, block=(self.nblocks,1,1), grid=(self.pixels_gpu.size//self.nblocks+1,1)) + self.context.pop() with self.lock: self.context.push() - t0 = time.time() - self.ray_trace_kernel(np.int32(self.pixels_gpu.size), self.origins_gpu, self.directions_gpu, self.pixels_gpu, block=(self.nblocks,1,1), grid=(self.pixels_gpu.size//self.nblocks+1,1)) - self.context.synchronize() - elapsed = time.time() - t0 - - print '\relapsed %f sec.' % elapsed, - sys.stdout.flush() - pygame.surfarray.blit_array(self.screen, self.pixels_gpu.get().reshape(self.size)) pygame.display.flip() self.context.pop() @@ -218,28 +233,23 @@ class Camera(Thread): clicked = False shift = False - current_layer = 0 + #current_layer = 0 while not done: self.clock.tick(20) if self.render and not clicked and not pygame.event.peek(KEYDOWN): - self.acquire_image() - self.render_image() + self.update() for event in pygame.event.get(): if event.type == MOUSEBUTTONDOWN: if event.button == 4: v = self.scale*np.cross(self.axis1,self.axis2)/10.0 self.translate(v) - self.clear_image() - self.update() if event.button == 5: v = -self.scale*np.cross(self.axis1,self.axis2)/10.0 self.translate(v) - self.clear_image() - self.update() if event.button == 1: clicked = True @@ -263,52 +273,41 @@ class Camera(Thread): if shift: v = mouse_direction*self.scale*length/float(self.width) self.translate(v) - self.clear_image() - self.update() else: phi = np.float32(2*np.pi*length/float(self.width)) n = rotate(mouse_direction, np.pi/2, -np.cross(self.axis1,self.axis2)) self.rotate(phi, n) - self.clear_image() - self.update() if event.type == KEYDOWN: if event.key == K_a: v = self.scale*self.axis1/10.0 self.translate(v) - self.clear_image() - self.update() if event.key == K_d: v = -self.scale*self.axis1/10.0 self.translate(v) - self.clear_image() - self.update() if event.key == K_w: v = self.scale*np.cross(self.axis1,self.axis2)/10.0 self.translate(v) - self.clear_image() - self.update() if event.key == K_s: v = -self.scale*np.cross(self.axis1,self.axis2)/10.0 self.translate(v) - self.clear_image() - self.update() if event.key == K_SPACE: v = self.scale*self.axis2/10.0 self.translate(v) - self.clear_image() - self.update() if event.key == K_LCTRL: v = -self.scale*self.axis2/10.0 self.translate(v) - self.clear_image() - self.update() + + if event.key == K_p: + self.clear_xyz_lookup() + self.update_xyz_lookup(self.point) + self.source_position = self.point if event.key == K_LSHIFT or event.key == K_RSHIFT: shift = True @@ -347,6 +346,9 @@ class Camera(Thread): self.screenshot() if event.key == K_F5: + if not hasattr(self, 'rng_states_gpu'): + self.initialize_render() + self.render = not self.render self.update() @@ -358,7 +360,6 @@ class Camera(Thread): done = True break - print pygame.display.quit() if __name__ == '__main__': diff --git a/color/Illuminanta.csv b/color/Illuminanta.csv new file mode 100644 index 0000000..3dba86c --- /dev/null +++ b/color/Illuminanta.csv @@ -0,0 +1,531 @@ +300,0.930483 +301,0.967643 +302,1.005970 +303,1.045490 +304,1.086230 +305,1.128210 +306,1.171470 +307,1.216020 +308,1.261880 +309,1.309100 +310,1.357690 +311,1.407680 +312,1.459100 +313,1.511980 +314,1.566330 +315,1.622190 +316,1.679590 +317,1.738550 +318,1.799100 +319,1.861270 +320,1.925080 +321,1.990570 +322,2.057760 +323,2.126670 +324,2.197340 +325,2.269800 +326,2.344060 +327,2.420170 +328,2.498140 +329,2.578010 +330,2.659810 +331,2.743550 +332,2.829280 +333,2.917010 +334,3.006780 +335,3.098610 +336,3.192530 +337,3.288570 +338,3.386760 +339,3.487120 +340,3.589680 +341,3.694470 +342,3.801520 +343,3.910850 +344,4.022500 +345,4.136480 +346,4.252820 +347,4.371560 +348,4.492720 +349,4.616310 +350,4.742380 +351,4.870950 +352,5.002040 +353,5.135680 +354,5.271890 +355,5.410700 +356,5.552130 +357,5.696220 +358,5.842980 +359,5.992440 +360,6.144620 +361,6.299550 +362,6.457240 +363,6.617740 +364,6.781050 +365,6.947200 +366,7.116210 +367,7.288110 +368,7.462920 +369,7.640660 +370,7.821350 +371,8.005010 +372,8.191670 +373,8.381340 +374,8.574040 +375,8.769800 +376,8.968640 +377,9.170560 +378,9.375610 +379,9.583780 +380,9.795100 +381,10.009600 +382,10.227300 +383,10.448100 +384,10.672200 +385,10.899600 +386,11.130200 +387,11.364000 +388,11.601200 +389,11.841600 +390,12.085300 +391,12.332400 +392,12.582800 +393,12.836600 +394,13.093800 +395,13.354300 +396,13.618200 +397,13.885500 +398,14.156300 +399,14.430400 +400,14.708000 +401,14.989100 +402,15.273600 +403,15.561600 +404,15.853000 +405,16.148000 +406,16.446400 +407,16.748400 +408,17.053800 +409,17.362800 +410,17.675300 +411,17.991300 +412,18.310800 +413,18.633900 +414,18.960500 +415,19.290700 +416,19.624400 +417,19.961700 +418,20.302600 +419,20.647000 +420,20.995000 +421,21.346500 +422,21.701600 +423,22.060300 +424,22.422500 +425,22.788300 +426,23.157700 +427,23.530700 +428,23.907200 +429,24.287300 +430,24.670900 +431,25.058100 +432,25.448900 +433,25.843200 +434,26.241100 +435,26.642500 +436,27.047500 +437,27.456000 +438,27.868100 +439,28.283600 +440,28.702700 +441,29.125300 +442,29.551500 +443,29.981100 +444,30.414200 +445,30.850800 +446,31.290900 +447,31.734500 +448,32.181500 +449,32.632000 +450,33.085900 +451,33.543200 +452,34.004000 +453,34.468200 +454,34.935800 +455,35.406800 +456,35.881100 +457,36.358800 +458,36.839900 +459,37.324300 +460,37.812100 +461,38.303100 +462,38.797500 +463,39.295100 +464,39.796000 +465,40.300200 +466,40.807600 +467,41.318200 +468,41.832000 +469,42.349100 +470,42.869300 +471,43.392600 +472,43.919200 +473,44.448800 +474,44.981600 +475,45.517400 +476,46.056300 +477,46.598300 +478,47.143300 +479,47.691300 +480,48.242300 +481,48.796300 +482,49.353300 +483,49.913200 +484,50.476000 +485,51.041800 +486,51.610400 +487,52.181800 +488,52.756100 +489,53.333200 +490,53.913200 +491,54.495800 +492,55.081300 +493,55.669400 +494,56.260300 +495,56.853900 +496,57.450100 +497,58.048900 +498,58.650400 +499,59.254500 +500,59.861100 +501,60.470300 +502,61.082000 +503,61.696200 +504,62.312800 +505,62.932000 +506,63.553500 +507,64.177500 +508,64.803800 +509,65.432500 +510,66.063500 +511,66.696800 +512,67.332400 +513,67.970200 +514,68.610200 +515,69.252500 +516,69.896900 +517,70.543500 +518,71.192200 +519,71.843000 +520,72.495900 +521,73.150800 +522,73.807700 +523,74.466600 +524,75.127500 +525,75.790300 +526,76.455100 +527,77.121700 +528,77.790200 +529,78.460500 +530,79.132600 +531,79.806500 +532,80.482100 +533,81.159500 +534,81.838600 +535,82.519300 +536,83.201700 +537,83.885600 +538,84.571200 +539,85.258400 +540,85.947000 +541,86.637200 +542,87.328800 +543,88.021900 +544,88.716500 +545,89.412400 +546,90.109700 +547,90.808300 +548,91.508200 +549,92.209500 +550,92.912000 +551,93.615700 +552,94.320600 +553,95.026700 +554,95.733900 +555,96.442300 +556,97.151800 +557,97.862300 +558,98.573900 +559,99.286400 +560,100.000000 +561,100.715000 +562,101.430000 +563,102.146000 +564,102.864000 +565,103.582000 +566,104.301000 +567,105.020000 +568,105.741000 +569,106.462000 +570,107.184000 +571,107.906000 +572,108.630000 +573,109.354000 +574,110.078000 +575,110.803000 +576,111.529000 +577,112.255000 +578,112.982000 +579,113.709000 +580,114.436000 +581,115.164000 +582,115.893000 +583,116.622000 +584,117.351000 +585,118.080000 +586,118.810000 +587,119.540000 +588,120.270000 +589,121.001000 +590,121.731000 +591,122.462000 +592,123.193000 +593,123.924000 +594,124.655000 +595,125.386000 +596,126.118000 +597,126.849000 +598,127.580000 +599,128.312000 +600,129.043000 +601,129.774000 +602,130.505000 +603,131.236000 +604,131.966000 +605,132.697000 +606,133.427000 +607,134.157000 +608,134.887000 +609,135.617000 +610,136.346000 +611,137.075000 +612,137.804000 +613,138.532000 +614,139.260000 +615,139.988000 +616,140.715000 +617,141.441000 +618,142.167000 +619,142.893000 +620,143.618000 +621,144.343000 +622,145.067000 +623,145.790000 +624,146.513000 +625,147.235000 +626,147.957000 +627,148.678000 +628,149.398000 +629,150.117000 +630,150.836000 +631,151.554000 +632,152.271000 +633,152.988000 +634,153.704000 +635,154.418000 +636,155.132000 +637,155.845000 +638,156.558000 +639,157.269000 +640,157.979000 +641,158.689000 +642,159.397000 +643,160.104000 +644,160.811000 +645,161.516000 +646,162.221000 +647,162.924000 +648,163.626000 +649,164.327000 +650,165.028000 +651,165.726000 +652,166.424000 +653,167.121000 +654,167.816000 +655,168.510000 +656,169.203000 +657,169.895000 +658,170.586000 +659,171.275000 +660,171.963000 +661,172.650000 +662,173.335000 +663,174.019000 +664,174.702000 +665,175.383000 +666,176.063000 +667,176.741000 +668,177.419000 +669,178.094000 +670,178.769000 +671,179.441000 +672,180.113000 +673,180.783000 +674,181.451000 +675,182.118000 +676,182.783000 +677,183.447000 +678,184.109000 +679,184.770000 +680,185.429000 +681,186.087000 +682,186.743000 +683,187.397000 +684,188.050000 +685,188.701000 +686,189.350000 +687,189.998000 +688,190.644000 +689,191.288000 +690,191.931000 +691,192.572000 +692,193.211000 +693,193.849000 +694,194.484000 +695,195.118000 +696,195.750000 +697,196.381000 +698,197.009000 +699,197.636000 +700,198.261000 +701,198.884000 +702,199.506000 +703,200.125000 +704,200.743000 +705,201.359000 +706,201.972000 +707,202.584000 +708,203.195000 +709,203.803000 +710,204.409000 +711,205.013000 +712,205.616000 +713,206.216000 +714,206.815000 +715,207.411000 +716,208.006000 +717,208.599000 +718,209.189000 +719,209.778000 +720,210.365000 +721,210.949000 +722,211.532000 +723,212.112000 +724,212.691000 +725,213.268000 +726,213.842000 +727,214.415000 +728,214.985000 +729,215.553000 +730,216.120000 +731,216.684000 +732,217.246000 +733,217.806000 +734,218.364000 +735,218.920000 +736,219.473000 +737,220.025000 +738,220.574000 +739,221.122000 +740,221.667000 +741,222.210000 +742,222.751000 +743,223.290000 +744,223.826000 +745,224.361000 +746,224.893000 +747,225.423000 +748,225.951000 +749,226.477000 +750,227.000000 +751,227.522000 +752,228.041000 +753,228.558000 +754,229.073000 +755,229.585000 +756,230.096000 +757,230.604000 +758,231.110000 +759,231.614000 +760,232.115000 +761,232.615000 +762,233.112000 +763,233.606000 +764,234.099000 +765,234.589000 +766,235.078000 +767,235.564000 +768,236.047000 +769,236.529000 +770,237.008000 +771,237.485000 +772,237.959000 +773,238.432000 +774,238.902000 +775,239.370000 +776,239.836000 +777,240.299000 +778,240.760000 +779,241.219000 +780,241.675000 +781,242.130000 +782,242.582000 +783,243.031000 +784,243.479000 +785,243.924000 +786,244.367000 +787,244.808000 +788,245.246000 +789,245.682000 +790,246.116000 +791,246.548000 +792,246.977000 +793,247.404000 +794,247.829000 +795,248.251000 +796,248.671000 +797,249.089000 +798,249.505000 +799,249.918000 +800,250.329000 +801,250.738000 +802,251.144000 +803,251.548000 +804,251.950000 +805,252.350000 +806,252.747000 +807,253.142000 +808,253.535000 +809,253.925000 +810,254.314000 +811,254.700000 +812,255.083000 +813,255.465000 +814,255.844000 +815,256.221000 +816,256.595000 +817,256.968000 +818,257.338000 +819,257.706000 +820,258.071000 +821,258.434000 +822,258.795000 +823,259.154000 +824,259.511000 +825,259.865000 +826,260.217000 +827,260.567000 +828,260.914000 +829,261.259000 +830,261.602000 \ No newline at end of file diff --git a/color/chromaticity.py b/color/chromaticity.py index fbd812a..859cc14 100644 --- a/color/chromaticity.py +++ b/color/chromaticity.py @@ -1,24 +1,14 @@ import numpy as np import os +import sys dir = os.path.split(os.path.realpath(__file__))[0] -f = open(dir + '/sbrgb10w.csv') +sys.path.append(dir + '/..') -color_map = [] -for line in f: - color_map.append([float(s) for s in line.split(',')]) +from tools import read_csv -f.close() - -color_map = np.array(color_map) - -# zero negative coefficients -color_map[color_map < 0] = 0 - -# normalize coefficients -for i in range(len(color_map)): - color_map[i,1:] /= np.sum(color_map[i,1:]) +color_map = read_csv(dir + '/ciexyz64_1.csv') def map_wavelength(wavelength): r = np.interp(wavelength, color_map[:,0], color_map[:,1]) diff --git a/color/scvle_1.csv b/color/scvle_1.csv new file mode 100644 index 0000000..0ab7def --- /dev/null +++ b/color/scvle_1.csv @@ -0,0 +1,401 @@ +380, 0.0005890000 +381, 0.0006650000 +382, 0.0007520000 +383, 0.0008540000 +384, 0.0009720000 +385, 0.0011080000 +386, 0.0012680000 +387, 0.0014530000 +388, 0.0016680000 +389, 0.0019180000 +390, 0.0022090000 +391, 0.0025470000 +392, 0.0029390000 +393, 0.0033940000 +394, 0.0039210000 +395, 0.0045300000 +396, 0.0052400000 +397, 0.0060500000 +398, 0.0069800000 +399, 0.0080600000 +400, 0.0092900000 +401, 0.0107000000 +402, 0.0123100000 +403, 0.0141300000 +404, 0.0161900000 +405, 0.0185200000 +406, 0.0211300000 +407, 0.0240500000 +408, 0.0273000000 +409, 0.0308900000 +410, 0.0348400000 +411, 0.0391600000 +412, 0.0439000000 +413, 0.0490000000 +414, 0.0545000000 +415, 0.0604000000 +416, 0.0668000000 +417, 0.0736000000 +418, 0.0808000000 +419, 0.0885000000 +420, 0.0966000000 +421, 0.1052000000 +422, 0.1141000000 +423, 0.1235000000 +424, 0.1334000000 +425, 0.1436000000 +426, 0.1541000000 +427, 0.1651000000 +428, 0.1764000000 +429, 0.1879000000 +430, 0.1998000000 +431, 0.2119000000 +432, 0.2243000000 +433, 0.2369000000 +434, 0.2496000000 +435, 0.2625000000 +436, 0.2755000000 +437, 0.2886000000 +438, 0.3017000000 +439, 0.3149000000 +440, 0.3281000000 +441, 0.3412000000 +442, 0.3543000000 +443, 0.3673000000 +444, 0.3803000000 +445, 0.3931000000 +446, 0.4060000000 +447, 0.4180000000 +448, 0.4310000000 +449, 0.4430000000 +450, 0.4550000000 +451, 0.4670000000 +452, 0.4790000000 +453, 0.4900000000 +454, 0.5020000000 +455, 0.5130000000 +456, 0.5240000000 +457, 0.5350000000 +458, 0.5460000000 +459, 0.5570000000 +460, 0.5670000000 +461, 0.5780000000 +462, 0.5880000000 +463, 0.5990000000 +464, 0.6100000000 +465, 0.6200000000 +466, 0.6310000000 +467, 0.6420000000 +468, 0.6530000000 +469, 0.6640000000 +470, 0.6760000000 +471, 0.6870000000 +472, 0.6990000000 +473, 0.7100000000 +474, 0.7220000000 +475, 0.7340000000 +476, 0.7450000000 +477, 0.7570000000 +478, 0.7690000000 +479, 0.7810000000 +480, 0.7930000000 +481, 0.8050000000 +482, 0.8170000000 +483, 0.8280000000 +484, 0.8400000000 +485, 0.8510000000 +486, 0.8620000000 +487, 0.8730000000 +488, 0.8840000000 +489, 0.8940000000 +490, 0.9040000000 +491, 0.9140000000 +492, 0.9230000000 +493, 0.9320000000 +494, 0.9410000000 +495, 0.9490000000 +496, 0.9570000000 +497, 0.9640000000 +498, 0.9700000000 +499, 0.9760000000 +500, 0.9820000000 +501, 0.9860000000 +502, 0.9900000000 +503, 0.9940000000 +504, 0.9970000000 +505, 0.9980000000 +506, 1.0000000000 +507, 1.0000000000 +508, 1.0000000000 +509, 0.9980000000 +510, 0.9970000000 +511, 0.9940000000 +512, 0.9900000000 +513, 0.9860000000 +514, 0.9810000000 +515, 0.9750000000 +516, 0.9680000000 +517, 0.9610000000 +518, 0.9530000000 +519, 0.9440000000 +520, 0.9350000000 +521, 0.9250000000 +522, 0.9150000000 +523, 0.9040000000 +524, 0.8920000000 +525, 0.8800000000 +526, 0.8670000000 +527, 0.8540000000 +528, 0.8400000000 +529, 0.8260000000 +530, 0.8110000000 +531, 0.7960000000 +532, 0.7810000000 +533, 0.7650000000 +534, 0.7490000000 +535, 0.7330000000 +536, 0.7170000000 +537, 0.7000000000 +538, 0.6830000000 +539, 0.6670000000 +540, 0.6500000000 +541, 0.6330000000 +542, 0.6160000000 +543, 0.5990000000 +544, 0.5810000000 +545, 0.5640000000 +546, 0.5480000000 +547, 0.5310000000 +548, 0.5140000000 +549, 0.4970000000 +550, 0.4810000000 +551, 0.4650000000 +552, 0.4480000000 +553, 0.4330000000 +554, 0.4170000000 +555, 0.4020000000 +556, 0.3864000000 +557, 0.3715000000 +558, 0.3569000000 +559, 0.3427000000 +560, 0.3288000000 +561, 0.3151000000 +562, 0.3018000000 +563, 0.2888000000 +564, 0.2762000000 +565, 0.2639000000 +566, 0.2519000000 +567, 0.2403000000 +568, 0.2291000000 +569, 0.2182000000 +570, 0.2076000000 +571, 0.1974000000 +572, 0.1876000000 +573, 0.1782000000 +574, 0.1690000000 +575, 0.1602000000 +576, 0.1517000000 +577, 0.1436000000 +578, 0.1358000000 +579, 0.1284000000 +580, 0.1212000000 +581, 0.1143000000 +582, 0.1078000000 +583, 0.1015000000 +584, 0.0956000000 +585, 0.0899000000 +586, 0.0845000000 +587, 0.0793000000 +588, 0.0745000000 +589, 0.0699000000 +590, 0.0655000000 +591, 0.0613000000 +592, 0.0574000000 +593, 0.0537000000 +594, 0.0502000000 +595, 0.0469000000 +596, 0.0438000000 +597, 0.0409000000 +598, 0.0381600000 +599, 0.0355800000 +600, 0.0331500000 +601, 0.0308700000 +602, 0.0287400000 +603, 0.0267400000 +604, 0.0248700000 +605, 0.0231200000 +606, 0.0214700000 +607, 0.0199400000 +608, 0.0185100000 +609, 0.0171800000 +610, 0.0159300000 +611, 0.0147700000 +612, 0.0136900000 +613, 0.0126900000 +614, 0.0117500000 +615, 0.0108800000 +616, 0.0100700000 +617, 0.0093200000 +618, 0.0086200000 +619, 0.0079700000 +620, 0.0073700000 +621, 0.0068200000 +622, 0.0063000000 +623, 0.0058200000 +624, 0.0053800000 +625, 0.0049700000 +626, 0.0045900000 +627, 0.0042400000 +628, 0.0039130000 +629, 0.0036130000 +630, 0.0033350000 +631, 0.0030790000 +632, 0.0028420000 +633, 0.0026230000 +634, 0.0024210000 +635, 0.0022350000 +636, 0.0020620000 +637, 0.0019030000 +638, 0.0017570000 +639, 0.0016210000 +640, 0.0014970000 +641, 0.0013820000 +642, 0.0012760000 +643, 0.0011780000 +644, 0.0010880000 +645, 0.0010050000 +646, 0.0009280000 +647, 0.0008570000 +648, 0.0007920000 +649, 0.0007320000 +650, 0.0006770000 +651, 0.0006260000 +652, 0.0005790000 +653, 0.0005360000 +654, 0.0004960000 +655, 0.0004590000 +656, 0.0004250000 +657, 0.0003935000 +658, 0.0003645000 +659, 0.0003377000 +660, 0.0003129000 +661, 0.0002901000 +662, 0.0002689000 +663, 0.0002493000 +664, 0.0002313000 +665, 0.0002146000 +666, 0.0001991000 +667, 0.0001848000 +668, 0.0001716000 +669, 0.0001593000 +670, 0.0001480000 +671, 0.0001375000 +672, 0.0001277000 +673, 0.0001187000 +674, 0.0001104000 +675, 0.0001026000 +676, 0.0000954000 +677, 0.0000888000 +678, 0.0000826000 +679, 0.0000769000 +680, 0.0000715000 +681, 0.0000666000 +682, 0.0000620000 +683, 0.0000578000 +684, 0.0000538000 +685, 0.0000501000 +686, 0.0000467000 +687, 0.0000436000 +688, 0.0000406000 +689, 0.0000378900 +690, 0.0000353300 +691, 0.0000329500 +692, 0.0000307500 +693, 0.0000287000 +694, 0.0000267900 +695, 0.0000250100 +696, 0.0000233600 +697, 0.0000218200 +698, 0.0000203800 +699, 0.0000190500 +700, 0.0000178000 +701, 0.0000166400 +702, 0.0000155600 +703, 0.0000145400 +704, 0.0000136000 +705, 0.0000127300 +706, 0.0000119100 +707, 0.0000111400 +708, 0.0000104300 +709, 0.0000097600 +710, 0.0000091400 +711, 0.0000085600 +712, 0.0000080200 +713, 0.0000075100 +714, 0.0000070400 +715, 0.0000066000 +716, 0.0000061800 +717, 0.0000058000 +718, 0.0000054400 +719, 0.0000051000 +720, 0.0000047800 +721, 0.0000044900 +722, 0.0000042100 +723, 0.0000039510 +724, 0.0000037090 +725, 0.0000034820 +726, 0.0000032700 +727, 0.0000030700 +728, 0.0000028840 +729, 0.0000027100 +730, 0.0000025460 +731, 0.0000023930 +732, 0.0000022500 +733, 0.0000021150 +734, 0.0000019890 +735, 0.0000018700 +736, 0.0000017590 +737, 0.0000016550 +738, 0.0000015570 +739, 0.0000014660 +740, 0.0000013790 +741, 0.0000012990 +742, 0.0000012230 +743, 0.0000011510 +744, 0.0000010840 +745, 0.0000010220 +746, 0.0000009620 +747, 0.0000009070 +748, 0.0000008550 +749, 0.0000008060 +750, 0.0000007600 +751, 0.0000007160 +752, 0.0000006750 +753, 0.0000006370 +754, 0.0000006010 +755, 0.0000005670 +756, 0.0000005350 +757, 0.0000005050 +758, 0.0000004770 +759, 0.0000004500 +760, 0.0000004250 +761, 0.0000004010 +762, 0.0000003790 +763, 0.0000003580 +764, 0.0000003382 +765, 0.0000003196 +766, 0.0000003021 +767, 0.0000002855 +768, 0.0000002699 +769, 0.0000002552 +770, 0.0000002413 +771, 0.0000002282 +772, 0.0000002159 +773, 0.0000002042 +774, 0.0000001932 +775, 0.0000001829 +776, 0.0000001731 +777, 0.0000001638 +778, 0.0000001551 +779, 0.0000001468 +780, 0.0000001390 \ No newline at end of file diff --git a/color/vl1924e_1.csv b/color/vl1924e_1.csv new file mode 100644 index 0000000..8b240c0 --- /dev/null +++ b/color/vl1924e_1.csv @@ -0,0 +1,471 @@ +360, 0.0000039170000 +361, 0.0000043935810 +362, 0.0000049296040 +363, 0.0000055321360 +364, 0.0000062082450 +365, 0.0000069650000 +366, 0.0000078132190 +367, 0.0000087673360 +368, 0.0000098398440 +369, 0.0000110432300 +370, 0.0000123900000 +371, 0.0000138864100 +372, 0.0000155572800 +373, 0.0000174429600 +374, 0.0000195837500 +375, 0.0000220200000 +376, 0.0000248396500 +377, 0.0000280412600 +378, 0.0000315310400 +379, 0.0000352152100 +380, 0.0000390000000 +381, 0.0000428264000 +382, 0.0000469146000 +383, 0.0000515896000 +384, 0.0000571764000 +385, 0.0000640000000 +386, 0.0000723442100 +387, 0.0000822122400 +388, 0.0000935081600 +389, 0.0001061361000 +390, 0.0001200000000 +391, 0.0001349840000 +392, 0.0001514920000 +393, 0.0001702080000 +394, 0.0001918160000 +395, 0.0002170000000 +396, 0.0002469067000 +397, 0.0002812400000 +398, 0.0003185200000 +399, 0.0003572667000 +400, 0.0003960000000 +401, 0.0004337147000 +402, 0.0004730240000 +403, 0.0005178760000 +404, 0.0005722187000 +405, 0.0006400000000 +406, 0.0007245600000 +407, 0.0008255000000 +408, 0.0009411600000 +409, 0.0010698800000 +410, 0.0012100000000 +411, 0.0013620910000 +412, 0.0015307520000 +413, 0.0017203680000 +414, 0.0019353230000 +415, 0.0021800000000 +416, 0.0024548000000 +417, 0.0027640000000 +418, 0.0031178000000 +419, 0.0035264000000 +420, 0.0040000000000 +421, 0.0045462400000 +422, 0.0051593200000 +423, 0.0058292800000 +424, 0.0065461600000 +425, 0.0073000000000 +426, 0.0080865070000 +427, 0.0089087200000 +428, 0.0097676800000 +429, 0.0106644300000 +430, 0.0116000000000 +431, 0.0125731700000 +432, 0.0135827200000 +433, 0.0146296800000 +434, 0.0157150900000 +435, 0.0168400000000 +436, 0.0180073600000 +437, 0.0192144800000 +438, 0.0204539200000 +439, 0.0217182400000 +440, 0.0230000000000 +441, 0.0242946100000 +442, 0.0256102400000 +443, 0.0269585700000 +444, 0.0283512500000 +445, 0.0298000000000 +446, 0.0313108300000 +447, 0.0328836800000 +448, 0.0345211200000 +449, 0.0362257100000 +450, 0.0380000000000 +451, 0.0398466700000 +452, 0.0417680000000 +453, 0.0437660000000 +454, 0.0458426700000 +455, 0.0480000000000 +456, 0.0502436800000 +457, 0.0525730400000 +458, 0.0549805600000 +459, 0.0574587200000 +460, 0.0600000000000 +461, 0.0626019700000 +462, 0.0652775200000 +463, 0.0680420800000 +464, 0.0709110900000 +465, 0.0739000000000 +466, 0.0770160000000 +467, 0.0802664000000 +468, 0.0836668000000 +469, 0.0872328000000 +470, 0.0909800000000 +471, 0.0949175500000 +472, 0.0990458400000 +473, 0.1033674000000 +474, 0.1078846000000 +475, 0.1126000000000 +476, 0.1175320000000 +477, 0.1226744000000 +478, 0.1279928000000 +479, 0.1334528000000 +480, 0.1390200000000 +481, 0.1446764000000 +482, 0.1504693000000 +483, 0.1564619000000 +484, 0.1627177000000 +485, 0.1693000000000 +486, 0.1762431000000 +487, 0.1835581000000 +488, 0.1912735000000 +489, 0.1994180000000 +490, 0.2080200000000 +491, 0.2171199000000 +492, 0.2267345000000 +493, 0.2368571000000 +494, 0.2474812000000 +495, 0.2586000000000 +496, 0.2701849000000 +497, 0.2822939000000 +498, 0.2950505000000 +499, 0.3085780000000 +500, 0.3230000000000 +501, 0.3384021000000 +502, 0.3546858000000 +503, 0.3716986000000 +504, 0.3892875000000 +505, 0.4073000000000 +506, 0.4256299000000 +507, 0.4443096000000 +508, 0.4633944000000 +509, 0.4829395000000 +510, 0.5030000000000 +511, 0.5235693000000 +512, 0.5445120000000 +513, 0.5656900000000 +514, 0.5869653000000 +515, 0.6082000000000 +516, 0.6293456000000 +517, 0.6503068000000 +518, 0.6708752000000 +519, 0.6908424000000 +520, 0.7100000000000 +521, 0.7281852000000 +522, 0.7454636000000 +523, 0.7619694000000 +524, 0.7778368000000 +525, 0.7932000000000 +526, 0.8081104000000 +527, 0.8224962000000 +528, 0.8363068000000 +529, 0.8494916000000 +530, 0.8620000000000 +531, 0.8738108000000 +532, 0.8849624000000 +533, 0.8954936000000 +534, 0.9054432000000 +535, 0.9148501000000 +536, 0.9237348000000 +537, 0.9320924000000 +538, 0.9399226000000 +539, 0.9472252000000 +540, 0.9540000000000 +541, 0.9602561000000 +542, 0.9660074000000 +543, 0.9712606000000 +544, 0.9760225000000 +545, 0.9803000000000 +546, 0.9840924000000 +547, 0.9874182000000 +548, 0.9903128000000 +549, 0.9928116000000 +550, 0.9949501000000 +551, 0.9967108000000 +552, 0.9980983000000 +553, 0.9991120000000 +554, 0.9997482000000 +555, 1.0000000000000 +556, 0.9998567000000 +557, 0.9993046000000 +558, 0.9983255000000 +559, 0.9968987000000 +560, 0.9950000000000 +561, 0.9926005000000 +562, 0.9897426000000 +563, 0.9864444000000 +564, 0.9827241000000 +565, 0.9786000000000 +566, 0.9740837000000 +567, 0.9691712000000 +568, 0.9638568000000 +569, 0.9581349000000 +570, 0.9520000000000 +571, 0.9454504000000 +572, 0.9384992000000 +573, 0.9311628000000 +574, 0.9234576000000 +575, 0.9154000000000 +576, 0.9070064000000 +577, 0.8982772000000 +578, 0.8892048000000 +579, 0.8797816000000 +580, 0.8700000000000 +581, 0.8598613000000 +582, 0.8493920000000 +583, 0.8386220000000 +584, 0.8275813000000 +585, 0.8163000000000 +586, 0.8047947000000 +587, 0.7930820000000 +588, 0.7811920000000 +589, 0.7691547000000 +590, 0.7570000000000 +591, 0.7447541000000 +592, 0.7324224000000 +593, 0.7200036000000 +594, 0.7074965000000 +595, 0.6949000000000 +596, 0.6822192000000 +597, 0.6694716000000 +598, 0.6566744000000 +599, 0.6438448000000 +600, 0.6310000000000 +601, 0.6181555000000 +602, 0.6053144000000 +603, 0.5924756000000 +604, 0.5796379000000 +605, 0.5668000000000 +606, 0.5539611000000 +607, 0.5411372000000 +608, 0.5283528000000 +609, 0.5156323000000 +610, 0.5030000000000 +611, 0.4904688000000 +612, 0.4780304000000 +613, 0.4656776000000 +614, 0.4534032000000 +615, 0.4412000000000 +616, 0.4290800000000 +617, 0.4170360000000 +618, 0.4050320000000 +619, 0.3930320000000 +620, 0.3810000000000 +621, 0.3689184000000 +622, 0.3568272000000 +623, 0.3447768000000 +624, 0.3328176000000 +625, 0.3210000000000 +626, 0.3093381000000 +627, 0.2978504000000 +628, 0.2865936000000 +629, 0.2756245000000 +630, 0.2650000000000 +631, 0.2547632000000 +632, 0.2448896000000 +633, 0.2353344000000 +634, 0.2260528000000 +635, 0.2170000000000 +636, 0.2081616000000 +637, 0.1995488000000 +638, 0.1911552000000 +639, 0.1829744000000 +640, 0.1750000000000 +641, 0.1672235000000 +642, 0.1596464000000 +643, 0.1522776000000 +644, 0.1451259000000 +645, 0.1382000000000 +646, 0.1315003000000 +647, 0.1250248000000 +648, 0.1187792000000 +649, 0.1127691000000 +650, 0.1070000000000 +651, 0.1014762000000 +652, 0.0961886400000 +653, 0.0911229600000 +654, 0.0862648500000 +655, 0.0816000000000 +656, 0.0771206400000 +657, 0.0728255200000 +658, 0.0687100800000 +659, 0.0647697600000 +660, 0.0610000000000 +661, 0.0573962100000 +662, 0.0539550400000 +663, 0.0506737600000 +664, 0.0475496500000 +665, 0.0445800000000 +666, 0.0417587200000 +667, 0.0390849600000 +668, 0.0365638400000 +669, 0.0342004800000 +670, 0.0320000000000 +671, 0.0299626100000 +672, 0.0280766400000 +673, 0.0263293600000 +674, 0.0247080500000 +675, 0.0232000000000 +676, 0.0218007700000 +677, 0.0205011200000 +678, 0.0192810800000 +679, 0.0181206900000 +680, 0.0170000000000 +681, 0.0159037900000 +682, 0.0148371800000 +683, 0.0138106800000 +684, 0.0128347800000 +685, 0.0119200000000 +686, 0.0110683100000 +687, 0.0102733900000 +688, 0.0095333110000 +689, 0.0088461570000 +690, 0.0082100000000 +691, 0.0076237810000 +692, 0.0070854240000 +693, 0.0065914760000 +694, 0.0061384850000 +695, 0.0057230000000 +696, 0.0053430590000 +697, 0.0049957960000 +698, 0.0046764040000 +699, 0.0043800750000 +700, 0.0041020000000 +701, 0.0038384530000 +702, 0.0035890990000 +703, 0.0033542190000 +704, 0.0031340930000 +705, 0.0029290000000 +706, 0.0027381390000 +707, 0.0025598760000 +708, 0.0023932440000 +709, 0.0022372750000 +710, 0.0020910000000 +711, 0.0019535870000 +712, 0.0018245800000 +713, 0.0017035800000 +714, 0.0015901870000 +715, 0.0014840000000 +716, 0.0013844960000 +717, 0.0012912680000 +718, 0.0012040920000 +719, 0.0011227440000 +720, 0.0010470000000 +721, 0.0009765896000 +722, 0.0009111088000 +723, 0.0008501332000 +724, 0.0007932384000 +725, 0.0007400000000 +726, 0.0006900827000 +727, 0.0006433100000 +728, 0.0005994960000 +729, 0.0005584547000 +730, 0.0005200000000 +731, 0.0004839136000 +732, 0.0004500528000 +733, 0.0004183452000 +734, 0.0003887184000 +735, 0.0003611000000 +736, 0.0003353835000 +737, 0.0003114404000 +738, 0.0002891656000 +739, 0.0002684539000 +740, 0.0002492000000 +741, 0.0002313019000 +742, 0.0002146856000 +743, 0.0001992884000 +744, 0.0001850475000 +745, 0.0001719000000 +746, 0.0001597781000 +747, 0.0001486044000 +748, 0.0001383016000 +749, 0.0001287925000 +750, 0.0001200000000 +751, 0.0001118595000 +752, 0.0001043224000 +753, 0.0000973356000 +754, 0.0000908458700 +755, 0.0000848000000 +756, 0.0000791466700 +757, 0.0000738580000 +758, 0.0000689160000 +759, 0.0000643026700 +760, 0.0000600000000 +761, 0.0000559818700 +762, 0.0000522256000 +763, 0.0000487184000 +764, 0.0000454474700 +765, 0.0000424000000 +766, 0.0000395610400 +767, 0.0000369151200 +768, 0.0000344486800 +769, 0.0000321481600 +770, 0.0000300000000 +771, 0.0000279912500 +772, 0.0000261135600 +773, 0.0000243602400 +774, 0.0000227246100 +775, 0.0000212000000 +776, 0.0000197785500 +777, 0.0000184528500 +778, 0.0000172168700 +779, 0.0000160645900 +780, 0.0000149900000 +781, 0.0000139872800 +782, 0.0000130515500 +783, 0.0000121781800 +784, 0.0000113625400 +785, 0.0000106000000 +786, 0.0000098858770 +787, 0.0000092173040 +788, 0.0000085923620 +789, 0.0000080091330 +790, 0.0000074657000 +791, 0.0000069595670 +792, 0.0000064879950 +793, 0.0000060486990 +794, 0.0000056393960 +795, 0.0000052578000 +796, 0.0000049017710 +797, 0.0000045697200 +798, 0.0000042601940 +799, 0.0000039717390 +800, 0.0000037029000 +801, 0.0000034521630 +802, 0.0000032183020 +803, 0.0000030003000 +804, 0.0000027971390 +805, 0.0000026078000 +806, 0.0000024312200 +807, 0.0000022665310 +808, 0.0000021130130 +809, 0.0000019699430 +810, 0.0000018366000 +811, 0.0000017122300 +812, 0.0000015962280 +813, 0.0000014880900 +814, 0.0000013873140 +815, 0.0000012934000 +816, 0.0000012058200 +817, 0.0000011241430 +818, 0.0000010480090 +819, 0.0000009770578 +820, 0.0000009109300 +821, 0.0000008492513 +822, 0.0000007917212 +823, 0.0000007380904 +824, 0.0000006881098 +825, 0.0000006415300 +826, 0.0000005980895 +827, 0.0000005575746 +828, 0.0000005198080 +829, 0.0000004846123 +830, 0.0000004518100 \ No newline at end of file diff --git a/detectors/__init__.py b/detectors/__init__.py index 199b4d7..5cffce4 100644 --- a/detectors/__init__.py +++ b/detectors/__init__.py @@ -24,6 +24,10 @@ def build_lbne_200kton(): def build_minilbne(): return build_lbne(radius/10, height/10, nstrings//10, pmts_per_string//10, endcap_spacing) +@buildable('microlbne') +def build_microlbne(): + return build_lbne(radius/40, height/40, nstrings//40, pmts_per_string//40, endcap_spacing) + @buildable('sno') def build_sno(): return build_sno_detector() diff --git a/geometry.py b/geometry.py index 1709a87..5fdfa1f 100644 --- a/geometry.py +++ b/geometry.py @@ -4,6 +4,7 @@ import pycuda.driver as cuda from pycuda import gpuarray from itertools import * from itertoolset import * +from camera import timeit # all material/surface properties are interpolated at these # wavelengths when they are sent to the gpu @@ -213,6 +214,7 @@ class Geometry(object): return len(self.solids)-1 + @timeit def build(self, bits=8, shift=3): offsets = [ (0,0) ] for solid in self.solids: @@ -352,18 +354,6 @@ class Geometry(object): if unique_bit_shifted_zvalues.size == 1: break - def load_colors(self, module): - self.colors = \ - np.concatenate([solid.color for solid in self.solids])[reorder] - - if not hasattr(self, 'colors_gpu'): - raise Exception('cannot load colors before geometry is loaded.') - - self.colors_gpu.set(self.colors.astype(np.uint32)) - - set_colors = module.get_function('set_colors') - set_colors(self.colors_gpu, block=(1,1,1), grid=(1,1)) - def load(self, module): """ Load the bounding volume hierarchy onto the GPU module, diff --git a/solids/__init__.py b/solids/__init__.py index 338725c..072c220 100644 --- a/solids/__init__.py +++ b/solids/__init__.py @@ -27,6 +27,12 @@ def build_12inch_pmt_with_lc(outer_material=water, theta=np.pi/8): pmt = build_12inch_pmt(outer_material, theta) return pmt + build_light_collector(pmt, a=lc_12inch_a, b=lc_12inch_b, d=lc_12inch_d, rmin=lc_12inch_rmin, rmax=lc_12inch_rmax) +@buildable('12inch_pmt_with_lc_hd') +def build_12inch_pmt_with_lc_hd(outer_material=water, theta=np.pi/32): + pmt = build_12inch_pmt(outer_material, theta) + return pmt + build_light_collector(pmt, a=lc_12inch_a, b=lc_12inch_b, d=lc_12inch_d, rmin=lc_12inch_rmin, rmax=lc_12inch_rmax, npoints=25) + + @buildable('8inch_pmt') def build_8inch_pmt(outer_material=water, theta=np.pi/8): return build_pmt(dir + '/sno_pmt.txt', 0.003, outer_material, theta) diff --git a/src/kernel.cu b/src/kernel.cu index a36ff91..be0087d 100644 --- a/src/kernel.cu +++ b/src/kernel.cu @@ -9,14 +9,12 @@ #include "photon.h" #include "alpha.h" -#define RED_WAVELENGTH 685 -#define BLUE_WAVELENGTH 445 -#define GREEN_WAVELENGTH 545 - __device__ void fAtomicAdd(float *addr, float data) { while (data) + { data = atomicExch(addr, data+atomicExch(addr, 0.0f)); + } } __device__ __noinline__ void to_diffuse(Photon &p, State &s, curandState &rng, const int &max_steps) @@ -28,9 +26,9 @@ __device__ __noinline__ void to_diffuse(Photon &p, State &s, curandState &rng, c int command; - command = fill_state(s, p, rng); + fill_state(s, p); - if (command == BREAK) + if (p.last_hit_triangle == -1) break; command = propagate_to_boundary(p, s, rng); @@ -87,104 +85,76 @@ __global__ void rotate(int nthreads, float3 *points, float phi, float3 axis) points[id] = rotate(points[id], phi, axis); } - __global__ void build_rgb_lookup(int nthreads, float3 *positions, float3 *directions, float3 *rgb_lookup1, float3 *rgb_lookup2, unsigned long long seed, int runs, int max_steps) +__global__ void update_xyz_lookup(int nthreads, int total_threads, int offset, float3 position, curandState *rng_states, float wavelength, float3 xyz, float3 *xyz_lookup1, float3 *xyz_lookup2, int max_steps) { - int id = blockIdx.x*blockDim.x + threadIdx.x; - - if (id >= nthreads) - return; - - curandState rng; - curand_init(seed, id, 0, &rng); - - Photon p0; - p0.position = positions[id]; - p0.direction = directions[id]; - p0.direction /= norm(p0.direction); - p0.polarization = uniform_sphere(&rng); - p0.time = 0.0f; - p0.history = 0x0; + int kernel_id = blockIdx.x*blockDim.x + threadIdx.x; + int id = kernel_id + offset; - float distance; - - int hit_triangle = intersect_mesh(p0.position, p0.direction, distance); - - if (hit_triangle != id) + if (kernel_id >= nthreads || id >= total_threads) return; - // note triangles from built geometry mesh are in order + curandState rng = rng_states[kernel_id]; - uint4 triangle_data = g_triangles[hit_triangle]; + uint4 triangle_data = g_triangles[id]; float3 v0 = g_vertices[triangle_data.x]; float3 v1 = g_vertices[triangle_data.y]; float3 v2 = g_vertices[triangle_data.z]; - float cos_theta = dot(normalize(cross(v1-v0, v2-v1)), -p0.direction); + float a = curand_uniform(&rng); + float b = uniform(&rng, 0.0f, (1.0f - a)); + float c = 1.0f - a - b; - if (cos_theta < 0.0f) - cos_theta = dot(-normalize(cross(v1-v0, v2-v1)), -p0.direction); + float3 direction = a*v0 + b*v1 + c*v2 - position; + direction /= norm(direction); - Photon p; - State s; + float distance; + int triangle_index = intersect_mesh(position, direction, distance); - for (int i=0; i < runs; i++) + if (triangle_index != id) { - p = p0; - p.wavelength = RED_WAVELENGTH; + rng_states[kernel_id] = rng; + return; + } - to_diffuse(p, s, rng, max_steps); + float cos_theta = dot(normalize(cross(v1-v0,v2-v1)),-direction); - if (p.history & REFLECT_DIFFUSE) - { - if (s.inside_to_outside) - { - fAtomicAdd(&rgb_lookup1[p.last_hit_triangle].x, cos_theta); - } - else - { - fAtomicAdd(&rgb_lookup2[p.last_hit_triangle].x, cos_theta); - } - } + if (cos_theta < 0.0f) + cos_theta = dot(-normalize(cross(v1-v0,v2-v1)),-direction); - p = p0; - p.wavelength = GREEN_WAVELENGTH; + Photon p; + p.position = position; + p.direction = direction; + p.wavelength = wavelength; + p.polarization = uniform_sphere(&rng); + p.last_hit_triangle = -1; + p.time = 0; + p.history = 0; - to_diffuse(p, s, rng, max_steps); + State s; + to_diffuse(p, s, rng, max_steps); - if (p.history & REFLECT_DIFFUSE) + if (p.history & REFLECT_DIFFUSE) + { + if (s.inside_to_outside) { - if (s.inside_to_outside) - { - fAtomicAdd(&rgb_lookup1[p.last_hit_triangle].y, cos_theta); - } - else - { - fAtomicAdd(&rgb_lookup2[p.last_hit_triangle].y, cos_theta); - } + fAtomicAdd(&xyz_lookup1[p.last_hit_triangle].x, cos_theta*xyz.x); + fAtomicAdd(&xyz_lookup1[p.last_hit_triangle].y, cos_theta*xyz.y); + fAtomicAdd(&xyz_lookup1[p.last_hit_triangle].z, cos_theta*xyz.z); } - - p = p0; - p.wavelength = BLUE_WAVELENGTH; - - to_diffuse(p, s, rng, max_steps); - - if (p.history & REFLECT_DIFFUSE) + else { - if (s.inside_to_outside) - { - fAtomicAdd(&rgb_lookup1[p.last_hit_triangle].z, cos_theta); - } - else - { - fAtomicAdd(&rgb_lookup2[p.last_hit_triangle].z, cos_theta); - } + fAtomicAdd(&xyz_lookup2[p.last_hit_triangle].x, cos_theta*xyz.x); + fAtomicAdd(&xyz_lookup2[p.last_hit_triangle].y, cos_theta*xyz.y); + fAtomicAdd(&xyz_lookup2[p.last_hit_triangle].z, cos_theta*xyz.z); } } -} // build_rgb_lookup + rng_states[kernel_id] = rng; -__global__ void render(int nthreads, curandState *rng_states, float3 *positions, float3 *directions, float *red, float *green, float *blue, float3 *rgb_lookup1, float3 *rgb_lookup2, int max_steps) +} // update_xyz_lookup + +__global__ void update_xyz_image(int nthreads, curandState *rng_states, float3 *positions, float3 *directions, float wavelength, float3 xyz, float3 *xyz_lookup1, float3 *xyz_lookup2, float3 *image, int nlookup_calls, int max_steps) { int id = blockIdx.x*blockDim.x + threadIdx.x; @@ -193,81 +163,69 @@ __global__ void render(int nthreads, curandState *rng_states, float3 *positions, curandState rng = rng_states[id]; - Photon p0; - p0.position = positions[id]; - p0.direction = directions[id]; - p0.direction /= norm(p0.direction); - p0.polarization = uniform_sphere(&rng); - p0.time = 0.0f; - p0.history = 0x0; + Photon p; + p.position = positions[id]; + p.direction = directions[id]; + p.direction /= norm(p.direction); + p.wavelength = wavelength; + p.polarization = uniform_sphere(&rng); + p.last_hit_triangle = -1; + p.time = 0; + p.history = 0; State s; - Photon p = p0; - p.wavelength = RED_WAVELENGTH; - to_diffuse(p, s, rng, max_steps); if (p.history & REFLECT_DIFFUSE) { if (s.inside_to_outside) { - red[id] = rgb_lookup1[p.last_hit_triangle].x; - } - else - { - red[id] = rgb_lookup2[p.last_hit_triangle].x; - } - } - else - { - red[id] = 0.0; - } - - p = p0; - p.wavelength = GREEN_WAVELENGTH; - - to_diffuse(p, s, rng, max_steps); - - if (p.history & REFLECT_DIFFUSE) - { - if (s.inside_to_outside) - { - green[id] = rgb_lookup1[p.last_hit_triangle].y; + image[id].x += xyz.x*xyz_lookup1[p.last_hit_triangle].x/nlookup_calls; + image[id].y += xyz.y*xyz_lookup1[p.last_hit_triangle].y/nlookup_calls; + image[id].z += xyz.z*xyz_lookup1[p.last_hit_triangle].z/nlookup_calls; } else { - green[id] = rgb_lookup2[p.last_hit_triangle].y; + image[id].x += xyz.x*xyz_lookup2[p.last_hit_triangle].x/nlookup_calls; + image[id].y += xyz.y*xyz_lookup2[p.last_hit_triangle].y/nlookup_calls; + image[id].z += xyz.z*xyz_lookup2[p.last_hit_triangle].z/nlookup_calls; } } - else - { - green[id] = 0.0; - } - - p = p0; - p.wavelength = BLUE_WAVELENGTH; - - to_diffuse(p, s, rng, max_steps); - - if (p.history & REFLECT_DIFFUSE) - { - if (s.inside_to_outside) - { - blue[id] = rgb_lookup1[p.last_hit_triangle].z; - } - else - { - blue[id] = rgb_lookup2[p.last_hit_triangle].z; - } - } - else - { - blue[id] = 0.0; - } - + rng_states[id] = rng; -} // render +} // update_xyz_image + +__global__ void process_image(int nthreads, float3 *image, int *pixels, int nimages) +{ + int id = blockIdx.x*blockDim.x + threadIdx.x; + + if (id >= nthreads) + return; + + float3 rgb = image[id]/nimages; + + if (rgb.x < 0.0f) + rgb.x = 0.0f; + if (rgb.y < 0.0f) + rgb.y = 0.0f; + if (rgb.z < 0.0f) + rgb.z = 0.0f; + + if (rgb.x > 1.0f) + rgb.x = 1.0f; + if (rgb.y > 1.0f) + rgb.y = 1.0f; + if (rgb.z > 1.0f) + rgb.z = 1.0f; + + unsigned int r = floorf(rgb.x*255.0f); + unsigned int g = floorf(rgb.y*255.0f); + unsigned int b = floorf(rgb.z*255.0f); + + pixels[id] = r << 16 | g << 8 | b; + +} // process_image /* Trace the rays starting at `positions` traveling in the direction `directions` to their intersection with the global mesh. If the ray @@ -322,9 +280,9 @@ __global__ void propagate(int nthreads, curandState *rng_states, float3 *positio int command; - command = fill_state(s, p, rng); + fill_state(s, p); - if (command == BREAK) + if (p.last_hit_triangle == -1) break; command = propagate_to_boundary(p, s, rng); diff --git a/src/photon.h b/src/photon.h index ad4c26c..f471866 100644 --- a/src/photon.h +++ b/src/photon.h @@ -19,8 +19,6 @@ struct Photon unsigned int history; int last_hit_triangle; - - //curandState rng; }; struct State @@ -51,14 +49,14 @@ enum enum {BREAK, CONTINUE, PASS}; // return value from propagate_to_boundary -__device__ int fill_state(State &s, Photon &p, curandState &rng) +__device__ void fill_state(State &s, Photon &p) { p.last_hit_triangle = intersect_mesh(p.position, p.direction, s.distance_to_boundary, p.last_hit_triangle); if (p.last_hit_triangle == -1) { p.history |= NO_HIT; - return BREAK; + return; } uint4 triangle_data = g_triangles[p.last_hit_triangle]; @@ -98,8 +96,6 @@ __device__ int fill_state(State &s, Photon &p, curandState &rng) s.absorption_length = interp_property(p.wavelength, material1.absorption_length); s.scattering_length = interp_property(p.wavelength, material1.scattering_length); - return PASS; - } // fill_state __device__ void rayleigh_scatter(Photon &p, curandState &rng) -- cgit From 2b1a23e20c08423dd6f0d290ec87bef2f836fbe1 Mon Sep 17 00:00:00 2001 From: Anthony LaTorre Date: Sat, 30 Jul 2011 21:03:45 -0400 Subject: you can rotate just the camera by holding the control key and take movies by pressing the m key. --- camera.py | 94 +++++++++++++++++++++++++++++++++++++++++++++++------- solids/__init__.py | 4 +-- src/kernel.cu | 24 ++++++++++---- 3 files changed, 102 insertions(+), 20 deletions(-) diff --git a/camera.py b/camera.py index 148ddf9..a39060d 100644 --- a/camera.py +++ b/camera.py @@ -16,6 +16,10 @@ from pycuda import gpuarray from pycuda.characterize import sizeof import pycuda.driver as cuda +from subprocess import call +import shutil +import tempfile + def timeit(func): def f(*args, **kwargs): t0 = time.time() @@ -24,6 +28,20 @@ def timeit(func): print '%s elapsed in %s().' % (datetime.timedelta(seconds=elapsed), func.__name__) return f +def encode_movie(dir): + root, ext = 'movie', 'avi' + for i in count(): + filename = '.'.join([root + str(i).zfill(5), ext]) + + if not os.path.exists(filename): + break + + call(['mencoder', 'mf://' + dir + '/*.png', '-mf', 'fps=10', '-o', filename, '-ovc', 'xvid', '-xvidencopts', 'bitrate=3000']) + + shutil.rmtree(dir) + + print 'movie saved to %s.' % filename + def get_rays(position, size = (800, 600), film_size = (0.035, 0.024), focal_length=0.05): """ Generate ray positions and directions from a pinhole camera facing the negative y direction. @@ -72,6 +90,7 @@ class Camera(Thread): self.context.push() self.ray_trace_kernel = self.module.get_function('ray_trace') self.rotate_kernel = self.module.get_function('rotate') + self.rotate_around_point_kernel = self.module.get_function('rotate_around_point') self.translate_kernel = self.module.get_function('translate') self.update_xyz_lookup_kernel = self.module.get_function('update_xyz_lookup') self.update_xyz_image_kernel = self.module.get_function('update_xyz_image') @@ -97,6 +116,9 @@ class Camera(Thread): self.pixels_gpu = gpuarray.zeros(self.width*self.height, dtype=np.int32) self.context.pop() + self.movie = False + self.movie_index = 0 + self.movie_dir = None self.render = False @timeit @@ -167,11 +189,11 @@ class Camera(Thread): self.process_image_kernel(np.int32(self.width*self.height), self.image_gpu, self.pixels_gpu, np.int32(self.nimages), block=(self.nblocks,1,1), grid=((self.width*self.height)//self.nblocks+1,1)) self.context.pop() - def screenshot(self, dir=''): + def screenshot(self, dir='', start=0): root, ext = 'screenshot', 'png' - for i in count(): - filename = os.path.join(dir, '.'.join([root + str(i).zfill(4), ext])) + for i in count(start): + filename = os.path.join(dir, '.'.join([root + str(i).zfill(5), ext])) if not os.path.exists(filename): break @@ -182,8 +204,8 @@ class Camera(Thread): def rotate(self, phi, n): with self.lock: self.context.push() - self.rotate_kernel(np.int32(self.pixels_gpu.size), self.origins_gpu, phi, gpuarray.vec.make_float3(*n), block=(self.nblocks,1,1), grid=(self.pixels_gpu.size//self.nblocks+1,1)) - self.rotate_kernel(np.int32(self.pixels_gpu.size), self.directions_gpu, phi, gpuarray.vec.make_float3(*n), block=(self.nblocks,1,1), grid=(self.pixels_gpu.size//self.nblocks+1,1)) + self.rotate_kernel(np.int32(self.pixels_gpu.size), self.origins_gpu, np.float32(phi), gpuarray.vec.make_float3(*n), block=(self.nblocks,1,1), grid=(self.pixels_gpu.size//self.nblocks+1,1)) + self.rotate_kernel(np.int32(self.pixels_gpu.size), self.directions_gpu, np.float32(phi), gpuarray.vec.make_float3(*n), block=(self.nblocks,1,1), grid=(self.pixels_gpu.size//self.nblocks+1,1)) self.point = rotate(self.point, phi, n) self.axis1 = rotate(self.axis1, phi, n) @@ -195,6 +217,21 @@ class Camera(Thread): self.update() + def rotate_around_point(self, phi, n, point): + 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)) + self.rotate_kernel(np.int32(self.directions_gpu.size), self.directions_gpu, np.float32(phi), gpuarray.vec.make_float3(*n), block=(self.nblocks,1,1), grid=(self.directions_gpu.size//self.nblocks+1,1)) + self.context.pop() + + self.axis1 = rotate(self.axis1, phi, n) + self.axis2 = rotate(self.axis2, phi, n) + + if self.render: + self.clear_image() + + self.update() + def translate(self, v): with self.lock: self.context.push() @@ -210,7 +247,7 @@ class Camera(Thread): def update(self): if self.render: - while self.nlookup_calls < 10: + while self.nlookup_calls < 100: self.update_xyz_lookup(self.source_position) self.update_image() self.process_image() @@ -226,12 +263,17 @@ class Camera(Thread): pygame.display.flip() self.context.pop() + if self.movie: + self.screenshot(self.movie_dir, self.movie_index) + self.movie_index += 1 + def run(self): self.update() done = False clicked = False shift = False + ctrl = False #current_layer = 0 @@ -277,7 +319,10 @@ class Camera(Thread): phi = np.float32(2*np.pi*length/float(self.width)) n = rotate(mouse_direction, np.pi/2, -np.cross(self.axis1,self.axis2)) - self.rotate(phi, n) + if ctrl: + self.rotate_around_point(phi, n, self.point) + else: + self.rotate(phi, n) if event.type == KEYDOWN: if event.key == K_a: @@ -300,18 +345,26 @@ class Camera(Thread): v = self.scale*self.axis2/10.0 self.translate(v) - if event.key == K_LCTRL: - v = -self.scale*self.axis2/10.0 - self.translate(v) + # if event.key == K_LCTRL: + # v = -self.scale*self.axis2/10.0 + # self.translate(v) - if event.key == K_p: + if event.key == K_F6: self.clear_xyz_lookup() - self.update_xyz_lookup(self.point) + self.clear_image() + self.source_position = self.point + + if event.key == K_p: + for i in range(100): + self.update_xyz_lookup(self.point) self.source_position = self.point if event.key == K_LSHIFT or event.key == K_RSHIFT: shift = True + if event.key == K_LCTRL or event.key == K_RCTRL: + ctrl = True + if event.key == K_ESCAPE: done = True break @@ -350,16 +403,33 @@ class Camera(Thread): self.initialize_render() self.render = not self.render + self.clear_image() self.update() + if event.key == K_m: + if self.movie: + encode_movie(self.movie_dir) + self.movie_dir = None + self.movie = False + else: + self.movie_index = 0 + self.movie_dir = tempfile.mkdtemp() + self.movie = True + if event.type == KEYUP: if event.key == K_LSHIFT or event.key == K_RSHIFT: shift = False + if event.key == K_LCTRL or event.key == K_RCTRL: + ctrl = False + if event.type == pygame.QUIT: done = True break + if self.movie: + encode_movie(self.movie_dir) + pygame.display.quit() if __name__ == '__main__': diff --git a/solids/__init__.py b/solids/__init__.py index 072c220..d37af7e 100644 --- a/solids/__init__.py +++ b/solids/__init__.py @@ -28,9 +28,9 @@ def build_12inch_pmt_with_lc(outer_material=water, theta=np.pi/8): return pmt + build_light_collector(pmt, a=lc_12inch_a, b=lc_12inch_b, d=lc_12inch_d, rmin=lc_12inch_rmin, rmax=lc_12inch_rmax) @buildable('12inch_pmt_with_lc_hd') -def build_12inch_pmt_with_lc_hd(outer_material=water, theta=np.pi/32): +def build_12inch_pmt_with_lc_hd(outer_material=water, theta=np.pi/64): pmt = build_12inch_pmt(outer_material, theta) - return pmt + build_light_collector(pmt, a=lc_12inch_a, b=lc_12inch_b, d=lc_12inch_d, rmin=lc_12inch_rmin, rmax=lc_12inch_rmax, npoints=25) + return pmt + build_light_collector(pmt, a=lc_12inch_a, b=lc_12inch_b, d=lc_12inch_d, rmin=lc_12inch_rmin, rmax=lc_12inch_rmax, npoints=100) @buildable('8inch_pmt') diff --git a/src/kernel.cu b/src/kernel.cu index be0087d..82efe88 100644 --- a/src/kernel.cu +++ b/src/kernel.cu @@ -62,27 +62,39 @@ __device__ __noinline__ void to_diffuse(Photon &p, State &s, curandState &rng, c extern "C" { -/* Translate `points` by the vector `v` */ -__global__ void translate(int nthreads, float3 *points, float3 v) +/* Translate the points `a` by the vector `v` */ +__global__ void translate(int nthreads, float3 *a, float3 v) { int id = blockIdx.x*blockDim.x + threadIdx.x; if (id >= nthreads) return; - points[id] += v; + a[id] += v; } -/* Rotate `points` through an angle `phi` counter-clockwise about the +/* Rotate the points `a` through an angle `phi` counter-clockwise about the axis `axis` (when looking towards +infinity). */ -__global__ void rotate(int nthreads, float3 *points, float phi, float3 axis) +__global__ void rotate(int nthreads, float3 *a, float phi, float3 axis) { int id = blockIdx.x*blockDim.x + threadIdx.x; if (id >= nthreads) return; - points[id] = rotate(points[id], phi, axis); + a[id] = rotate(a[id], phi, axis); +} + +__global__ void rotate_around_point(int nthreads, float3 *a, float phi, float3 axis, float3 point) +{ + int id = blockIdx.x*blockDim.x + threadIdx.x; + + if (id >= nthreads) + return; + + a[id] -= point; + a[id] = rotate(a[id], phi, axis); + a[id] += point; } __global__ void update_xyz_lookup(int nthreads, int total_threads, int offset, float3 position, curandState *rng_states, float wavelength, float3 xyz, float3 *xyz_lookup1, float3 *xyz_lookup2, int max_steps) -- cgit From 42d2241948e04953788ce82bbeef25d0cba584fb Mon Sep 17 00:00:00 2001 From: Anthony LaTorre Date: Sat, 30 Jul 2011 21:06:46 -0400 Subject: reduce number of lookup calls. --- camera.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/camera.py b/camera.py index a39060d..37dc072 100644 --- a/camera.py +++ b/camera.py @@ -247,7 +247,7 @@ class Camera(Thread): 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() @@ -355,7 +355,7 @@ class Camera(Thread): self.source_position = self.point if event.key == K_p: - for i in range(100): + for i in range(10): self.update_xyz_lookup(self.point) self.source_position = self.point -- cgit