Main Page   Namespace List   Class Hierarchy   Alphabetical List   Compound List   File List   Namespace Members   Compound Members   File Members   Related Pages  

OptiXShaders.cu

Go to the documentation of this file.
00001 /***************************************************************************
00002  *cr
00003  *cr            (C) Copyright 1995-2019 The Board of Trustees of the
00004  *cr                        University of Illinois
00005  *cr                         All Rights Reserved
00006  *cr
00007  ***************************************************************************/
00008 
00009 /***************************************************************************
00010 * RCS INFORMATION:
00011 *
00012 *      $RCSfile: OptiXShaders.cu,v $
00013 *      $Author: johns $      $Locker:  $               $State: Exp $
00014 *      $Revision: 1.176 $         $Date: 2021/05/14 22:15:01 $
00015 *
00016 ***************************************************************************/
00080 #include <optix.h>
00081 #include <optixu/optixu_math_namespace.h>
00082 #include <optixu/optixu_matrix_namespace.h>
00083 #include <optixu/optixu_aabb_namespace.h>
00084 #include "OptiXShaders.h"
00085 
00086 // XXX special shader modifications used in the VMD renderings
00087 // included in the NSF CADENS "Birth of Planet Earth" movie
00088 // #define CADENSMOVIE 1
00089 
00090 // Runtime-based RT output coloring
00091 //#define ORT_TIME_COLORING 1
00092 // normalize runtime color by max acceptable pixel runtime 
00093 #define ORT_TIME_NORMALIZATION (1e-9 / 2.5f)
00094 // the default time coloring method averages contributions,
00095 // but we often want to know what the worst time was rather than
00096 // the average. 
00097 #define ORT_TIME_COMBINE_MAX 1
00098 
00099 // Macros related to ray origin epsilon stepping to prevent 
00100 // self-intersections with the surface we're leaving
00101 // This is a cheesy way of avoiding self-intersection 
00102 // but it ameliorates the problem.
00103 // At present without the ray origin stepping, even a simple
00104 // STMV test case will exhibit self-intersection artifacts that
00105 // are particularly obvious in shadowing for direct lighting.
00106 // Since changing the scene epsilon even to large values does not 
00107 // always cure the problem, this workaround is still required.
00108 #define ORT_USE_RAY_STEP       1
00109 #define ORT_TRANS_USE_INCIDENT 1
00110 #define ORT_RAY_STEP           N*scene_epsilon*4.0f
00111 #define ORT_RAY_STEP2          ray.direction*scene_epsilon*4.0f
00112 
00113 // reverse traversal of any-hit rays for shadows/AO
00114 //
00115 // XXX The macro definition of USE_REVERSE_SHADOW_RAYS is now located in 
00116 //     OptiXShaders.h, since both the shader and the host code need to
00117 //     know if it is defined...
00118 //
00119 // XXX We currently hard-code REVERSE_RAY_LENGTH in such a way that
00120 //     it works well for scenes that fall within the VMD view volume,
00121 //     given the relationship between the model and camera coordinate
00122 //     systems, but this would be best computed by the diagonal of the 
00123 //     AABB for the full scene, and then scaled into camera coordinates.
00124 //     The REVERSE_RAY_STEP size is computed to avoid self intersection 
00125 //     with the surface we're shading.
00126 #define REVERSE_RAY_STEP       (scene_epsilon*10.0f)
00127 #define REVERSE_RAY_LENGTH     3.0f
00128 
00129 
00130 // Macros to enable particular ray-geometry intersection variants that
00131 // optimize for speed, or some combination of speed and accuracy
00132 #define ORT_USE_SPHERES_HEARNBAKER 1
00133 
00134 
00135 using namespace optix;
00136 
00137 //
00138 // TEA, a tiny encryption algorithm.
00139 // D. Wheeler and R. Needham, 2nd Intl. Workshop Fast Software Encryption, 
00140 // LNCS, pp. 363-366, 1994.
00141 //
00142 // GPU Random Numbers via the Tiny Encryption Algorithm
00143 // F. Zafar, M. Olano, and A. Curtis.
00144 // HPG '10 Proceedings of the Conference on High Performance Graphics,
00145 // pp. 133-141, 2010.  
00146 //
00147 template<unsigned int N> static __host__ __device__ __inline__ 
00148 unsigned int tea(unsigned int val0, unsigned int val1) {
00149   unsigned int v0 = val0;
00150   unsigned int v1 = val1;
00151   unsigned int s0 = 0;
00152 
00153   for( unsigned int n = 0; n < N; n++ ) {
00154     s0 += 0x9e3779b9;
00155     v0 += ((v1<<4)+0xa341316c)^(v1+s0)^((v1>>5)+0xc8013ea4);
00156     v1 += ((v0<<4)+0xad90777d)^(v0+s0)^((v0>>5)+0x7e95761e);
00157   }
00158 
00159   return v0;
00160 }
00161 
00162 
00163 // output image framebuffer, accumulation buffer, and ray launch parameters
00164 rtBuffer<uchar4, 2> framebuffer;
00165 rtBuffer<float4, 2> accumulation_buffer;
00166 rtDeclareVariable(float, accumulation_normalization_factor, , );
00167 rtDeclareVariable(uint2, launch_index, rtLaunchIndex, );
00168 rtDeclareVariable(uint2, launch_dim, rtLaunchDim, );
00169 #if defined(VMDOPTIX_PROGRESSIVEAPI)
00170 rtDeclareVariable(unsigned int, progressiveSubframeIndex, rtSubframeIndex, );
00171 #endif
00172 rtDeclareVariable(unsigned int, accumCount, , );
00173 rtDeclareVariable(int, progressive_enabled, , );
00174 
00175 #if defined(ORT_RAYSTATS)
00176 // input/output ray statistics buffer
00177 rtBuffer<uint4, 2> raystats1_buffer; // x=prim, y=shad-dir, z=shad-ao, w=miss
00178 rtBuffer<uint4, 2> raystats2_buffer; // x=trans, y=trans-skip, z=?, w=refl
00179 #endif
00180 
00181 // epsilon value to use to avoid self-intersection 
00182 rtDeclareVariable(float, scene_epsilon, , );
00183 
00184 // max ray recursion depth
00185 rtDeclareVariable(unsigned int, max_depth, , );
00186 
00187 // max number of transparent surfaces (max_trans <= max_depth)
00188 rtDeclareVariable(unsigned int, max_trans, , );
00189 
00190 // XXX global interpolation coordinate for experimental  animated 
00191 // representations that loop over some fixed sequence of motion 
00192 // over the domain [0:1].
00193 rtDeclareVariable(float, anim_interp, , );
00194 
00195 // shadow rendering mode
00196 rtDeclareVariable(int, shadows_enabled, , );
00197 rtDeclareVariable(int, aa_samples, , );
00198 
00199 // ambient occlusion sample counts and scaling factors
00200 rtDeclareVariable(int, ao_samples, , );
00201 rtDeclareVariable(float, ao_ambient, , );
00202 rtDeclareVariable(float, ao_direct, , );
00203 rtDeclareVariable(float, ao_maxdist, , ); 
00204 
00205 // background color and/or background gradient
00206 rtDeclareVariable(float3, scene_bg_color, , );
00207 rtDeclareVariable(float3, scene_bg_color_grad_top, , );
00208 rtDeclareVariable(float3, scene_bg_color_grad_bot, , );
00209 rtDeclareVariable(float3, scene_gradient, , );
00210 rtDeclareVariable(float, scene_gradient_topval, , );
00211 rtDeclareVariable(float, scene_gradient_botval, , );
00212 rtDeclareVariable(float, scene_gradient_invrange, , );
00213 
00214 // VR HMD fade+clipping plane/sphere
00215 rtDeclareVariable(int, clipview_mode, , );
00216 rtDeclareVariable(float, clipview_start, , );
00217 rtDeclareVariable(float, clipview_end, , );
00218 
00219 // VR HMD headlight
00220 rtDeclareVariable(int, headlight_mode, , );
00221 
00222 // fog state
00223 rtDeclareVariable(int, fog_mode, , );
00224 rtDeclareVariable(float, fog_start, , );
00225 rtDeclareVariable(float, fog_end, , );
00226 rtDeclareVariable(float, fog_density, , );
00227 
00228 // camera parameters
00229 rtDeclareVariable(float,  cam_zoom, , );
00230 rtDeclareVariable(float3, cam_pos, , );
00231 rtDeclareVariable(float3, cam_U, , );
00232 rtDeclareVariable(float3, cam_V, , );
00233 rtDeclareVariable(float3, cam_W, , );
00234 
00235 // stereoscopic camera parameters
00236 rtDeclareVariable(float,  cam_stereo_eyesep, , );
00237 rtDeclareVariable(float,  cam_stereo_convergence_dist, , );
00238 
00239 // camera depth-of-field parameters
00240 rtDeclareVariable(float,  cam_dof_aperture_rad, , );
00241 rtDeclareVariable(float,  cam_dof_focal_dist, , );
00242 
00243 // various shading related per-ray state
00244 rtDeclareVariable(unsigned int, radiance_ray_type, , );
00245 rtDeclareVariable(unsigned int, shadow_ray_type, , );
00246 rtDeclareVariable(rtObject, root_object, , );
00247 rtDeclareVariable(rtObject, root_shadower, , );
00248 
00249 rtDeclareVariable(optix::Ray, ray, rtCurrentRay, );
00250 rtDeclareVariable(float, t_hit, rtIntersectionDistance, );
00251 //rtDeclareVariable(float3, texcoord, attribute texcoord, );
00252 
00253 //
00254 // Classic API
00255 //
00256 rtDeclareVariable(float3, geometric_normal, attribute geometric_normal, );
00257 rtDeclareVariable(float3, shading_normal, attribute shading_normal, );
00258 // object color assigned at intersection time
00259 rtDeclareVariable(float3, obj_color, attribute obj_color, );
00260 
00261 #if defined(ORT_USERTXAPIS)
00262 //
00263 // OptiX RTX hardware triangle API
00264 //
00265 rtDeclareVariable(float2, barycentrics, attribute rtTriangleBarycentrics, );
00266 rtBuffer<uint4> normalBuffer; // packed normals: ng [n0 n1 n2]
00267 rtBuffer<uchar4> colorBuffer; // unsigned char color representation
00268 rtDeclareVariable(int, has_vertex_normals, , );
00269 rtDeclareVariable(int, has_vertex_colors, , );
00270 #endif
00271 
00272 
00273 struct PerRayData_radiance {
00274   float3 result;         // final shaded surface color
00275   float alpha;           // alpha value to back-propagate to framebuffer
00276   float importance;      // importance of recursive ray tree
00277   unsigned int depth;    // current recursion depth
00278   unsigned int transcnt; // transmission ray surface count/depth
00279 };
00280 
00281 struct PerRayData_shadow {
00282   float3 attenuation;
00283 };
00284 
00285 rtDeclareVariable(PerRayData_radiance, prd_radiance, rtPayload, );
00286 rtDeclareVariable(PerRayData_radiance, prd, rtPayload, );
00287 rtDeclareVariable(PerRayData_shadow, prd_shadow, rtPayload, );
00288 
00289 // list of directional and positional lights
00290 #if defined(VMDOPTIX_LIGHTUSEROBJS)
00291 rtDeclareVariable(DirectionalLightList, dir_light_list, , );
00292 rtDeclareVariable(PositionalLightList, pos_light_list, , );
00293 #else
00294 rtBuffer<DirectionalLight> dir_lights;
00295 rtBuffer<DirectionalLight> pos_lights;
00296 #endif
00297 
00298 //
00299 // convert float3 rgb data to uchar4 with alpha channel set to 255.
00300 //
00301 static __device__ __inline__ uchar4 make_color_rgb4u(const float3& c) {
00302   return make_uchar4(static_cast<unsigned char>(__saturatef(c.x)*255.99f),  
00303                      static_cast<unsigned char>(__saturatef(c.y)*255.99f),  
00304                      static_cast<unsigned char>(__saturatef(c.z)*255.99f),  
00305                      255u);
00306 }
00307 
00308 //
00309 // convert float4 rgba data to uchar4
00310 //
00311 static __device__ __inline__ uchar4 make_color_rgb4u(const float4& c) {
00312   return make_uchar4(static_cast<unsigned char>(__saturatef(c.x)*255.99f),  
00313                      static_cast<unsigned char>(__saturatef(c.y)*255.99f),  
00314                      static_cast<unsigned char>(__saturatef(c.z)*255.99f),  
00315                      static_cast<unsigned char>(__saturatef(c.w)*255.99f));
00316 }
00317 
00318 
00319 //
00320 // OptiX miss programs for drawing the background color or 
00321 // background color gradient when no objects are hit
00322 //
00323 
00324 // Miss program for solid background
00325 RT_PROGRAM void miss_solid_bg() {
00326   // Fog overrides the background color if we're using
00327   // Tachyon radial fog, but not for OpenGL style fog.
00328   prd_radiance.result = scene_bg_color;
00329   prd_radiance.alpha = 0.0f; // alpha of background is 0.0f;
00330 #if defined(ORT_RAYSTATS)
00331   raystats1_buffer[launch_index].w++; // increment miss counter
00332 #endif
00333 }
00334 
00335 
00336 // Miss program for gradient background with perspective projection,
00337 // adapted from Tachyon
00338 // Fog overrides the background color if we're using
00339 // Tachyon radial fog, but not for OpenGL style fog.
00340 RT_PROGRAM void miss_gradient_bg_sky_sphere() {
00341   float IdotG = dot(ray.direction, scene_gradient);
00342   float val = (IdotG - scene_gradient_botval) * scene_gradient_invrange;
00343   val = __saturatef(val);
00344   float3 col = val * scene_bg_color_grad_top + 
00345                (1.0f - val) * scene_bg_color_grad_bot; 
00346   prd_radiance.result = col;
00347   prd_radiance.alpha = 0.0f; // alpha of background is 0.0f;
00348 #if defined(ORT_RAYSTATS)
00349   raystats1_buffer[launch_index].w++; // increment miss counter
00350 #endif
00351 }
00352 
00353 
00354 // Miss program for gradient background with orthographic projection,
00355 // adapted from Tachyon
00356 // Fog overrides the background color if we're using
00357 // Tachyon radial fog, but not for OpenGL style fog.
00358 RT_PROGRAM void miss_gradient_bg_sky_plane() {
00359   float IdotG = dot(ray.origin, scene_gradient);
00360   float val = (IdotG - scene_gradient_botval) * scene_gradient_invrange;
00361   val = __saturatef(val);
00362   float3 col = val * scene_bg_color_grad_top + 
00363                (1.0f - val) * scene_bg_color_grad_bot; 
00364   prd_radiance.result = col;
00365   prd_radiance.alpha = 0.0f; // alpha of background is 0.0f;
00366 #if defined(ORT_RAYSTATS)
00367   raystats1_buffer[launch_index].w++; // increment miss counter
00368 #endif
00369 }
00370 
00371 
00372 //
00373 // Various random number routines adapted from Tachyon
00374 //
00375 #define MYRT_RAND_MAX 4294967296.0f       /* Max random value from rt_rand  */
00376 #define MYRT_RAND_MAX_INV .00000000023283064365f   /* normalize rt_rand  */
00377 
00378 //
00379 // Quick and dirty 32-bit LCG random number generator [Fishman 1990]:
00380 //   A=1099087573 B=0 M=2^32
00381 //   Period: 10^9
00382 // Fastest gun in the west, but fails many tests after 10^6 samples,
00383 // and fails all statistics tests after 10^7 samples.
00384 // It fares better than the Numerical Recipes LCG.  This is the fastest
00385 // power of two rand, and has the best multiplier for 2^32, found by
00386 // brute force[Fishman 1990].  Test results:
00387 //   http://www.iro.umontreal.ca/~lecuyer/myftp/papers/testu01.pdf
00388 //   http://www.shadlen.org/ichbin/random/
00389 //
00390 static __device__ __inline__ 
00391 unsigned int myrt_rand(unsigned int &idum) {
00392   // on machines where int is 32-bits, no need to mask
00393   idum *= 1099087573;
00394   return idum;
00395 }
00396 
00397 
00398 // Generate an offset to jitter AA samples in the image plane, adapted
00399 // from the code in Tachyon
00400 static __device__ __inline__ 
00401 void jitter_offset2f(unsigned int &pval, float2 &xy) {
00402   xy.x = (myrt_rand(pval) * MYRT_RAND_MAX_INV) - 0.5f;
00403   xy.y = (myrt_rand(pval) * MYRT_RAND_MAX_INV) - 0.5f;
00404 }
00405 
00406 
00407 // Generate an offset to jitter DoF samples in the Circle of Confusion,
00408 // adapted from the code in Tachyon
00409 static __device__ __inline__ 
00410 void jitter_disc2f(unsigned int &pval, float2 &xy, float radius) {
00411 #if 1
00412   // Since the GPU RT currently uses super cheap/sleazy LCG RNGs,
00413   // it is best to avoid using sample picking, which can fail if
00414   // we use a multiply-only RNG and we hit a zero in the PRN sequence.
00415   // The special functions are slow, but have bounded runtime and 
00416   // minimal branch divergence.
00417   float   r=(myrt_rand(pval) * MYRT_RAND_MAX_INV);
00418   float phi=(myrt_rand(pval) * MYRT_RAND_MAX_INV) * 2.0f * M_PIf;
00419   __sincosf(phi, &xy.x, &xy.y); // fast approximation
00420   xy *= sqrtf(r) * radius;
00421 #else
00422   // Pick uniform samples that fall within the disc --
00423   // this scheme can hang in an endless loop if a poor quality
00424   // RNG is used and it gets stuck in a short PRN sub-sequence
00425   do { 
00426     xy.x = 2.0f * (myrt_rand(pval) * MYRT_RAND_MAX_INV) - 1.0f;
00427     xy.y = 2.0f * (myrt_rand(pval) * MYRT_RAND_MAX_INV) - 1.0f;
00428   } while ((xy.x*xy.x + xy.y*xy.y) > 1.0f);
00429   xy *= radius; 
00430 #endif
00431 }
00432 
00433 
00434 // Generate a randomly oriented ray, based on Tachyon implementation
00435 static __device__ __inline__ 
00436 void jitter_sphere3f(unsigned int &pval, float3 &dir) {
00437 #if 1
00438   //
00439   // Use GPU fast/approximate math routines
00440   //
00441   /* Archimedes' cylindrical projection scheme       */
00442   /* generate a point on a unit cylinder and project */
00443   /* back onto the sphere.  This approach is likely  */
00444   /* faster for SIMD hardware, despite the use of    */
00445   /* transcendental functions.                       */
00446   float u1 = myrt_rand(pval) * MYRT_RAND_MAX_INV;
00447   dir.z = 2.0f * u1 - 1.0f;
00448   float R = __fsqrt_rn(1.0f - dir.z*dir.z);  // fast approximation
00449   float u2 = myrt_rand(pval) * MYRT_RAND_MAX_INV;
00450   float phi = 2.0f * M_PIf * u2;
00451   float sinphi, cosphi;
00452   __sincosf(phi, &sinphi, &cosphi); // fast approximation
00453   dir.x = R * cosphi;
00454   dir.y = R * sinphi;
00455 #elif 1
00456   /* Archimedes' cylindrical projection scheme       */
00457   /* generate a point on a unit cylinder and project */
00458   /* back onto the sphere.  This approach is likely  */
00459   /* faster for SIMD hardware, despite the use of    */
00460   /* transcendental functions.                       */
00461   float u1 = myrt_rand(pval) * MYRT_RAND_MAX_INV;
00462   dir.z = 2.0f * u1 - 1.0f;
00463   float R = sqrtf(1.0f - dir.z*dir.z); 
00464 
00465   float u2 = myrt_rand(pval) * MYRT_RAND_MAX_INV;
00466   float phi = 2.0f * M_PIf * u2;
00467   float sinphi, cosphi;
00468   sincosf(phi, &sinphi, &cosphi);
00469   dir.x = R * cosphi;
00470   dir.y = R * sinphi;
00471 #else
00472   /* Marsaglia's uniform sphere sampling scheme           */
00473   /* In order to correctly sample a sphere, using rays    */
00474   /* generated randomly within a cube we must throw out   */
00475   /* direction vectors longer than 1.0, otherwise we'll   */
00476   /* oversample the corners of the cube relative to       */
00477   /* a true sphere.                                       */
00478   float len;
00479   float3 d;
00480   do {
00481     d.x = (myrt_rand(pval) * MYRT_RAND_MAX_INV) - 0.5f;
00482     d.y = (myrt_rand(pval) * MYRT_RAND_MAX_INV) - 0.5f;
00483     d.z = (myrt_rand(pval) * MYRT_RAND_MAX_INV) - 0.5f;
00484     len = dot(d, d);
00485   } while (len > 0.250f);
00486   float invlen = rsqrtf(len);
00487 
00488   /* finish normalizing the direction vector */
00489   dir = d * invlen;
00490 #endif
00491 }
00492 
00493 
00494 //
00495 // Device functions for clipping rays by geometric primitives
00496 // 
00497 
00498 // fade_start: onset of fading 
00499 //   fade_end: fully transparent, begin clipping of geometry
00500 __device__ void sphere_fade_and_clip(const float3 &hit_point, 
00501                                      const float3 &cam_pos,
00502                                      float fade_start, float fade_end,
00503                                      float &alpha) {
00504   float camdist = length(hit_point - cam_pos);
00505 
00506   // we can omit the distance test since alpha modulation value is clamped
00507   // if (1 || camdist < fade_start) {
00508     float fade_len = fade_start - fade_end;
00509     alpha *= __saturatef((camdist - fade_start) / fade_len);
00510   // }
00511 }
00512 
00513 
00514 __device__ void ray_sphere_clip_interval(const optix::Ray &ray, float3 center,
00515                                          float rad, float2 &tinterval) {
00516   float3 V = center - ray.origin;
00517   float b = dot(V, ray.direction);
00518   float disc = b*b + rad*rad - dot(V, V);
00519 
00520   // if the discriminant is positive, the ray hits...
00521   if (disc > 0.0f) {
00522     disc = sqrtf(disc);
00523     tinterval.x = b-disc;
00524     tinterval.y = b+disc;
00525   } else {
00526     tinterval.x = -RT_DEFAULT_MAX; 
00527     tinterval.y =  RT_DEFAULT_MAX; 
00528   }
00529 }
00530 
00531 
00532 __device__ void clip_ray_by_plane(optix::Ray &ray, const float4 plane) {
00533   float3 n = make_float3(plane);                                              
00534   float dt = dot(ray.direction, n);                                            
00535   float t = (-plane.w - dot(n, ray.origin))/dt;                                 
00536   if(t > ray.tmin && t < ray.tmax) {                                          
00537     if (dt <= 0) {                                                              
00538       ray.tmax = t;                                                             
00539     } else {                                                                    
00540       ray.tmin = t;                                                             
00541     }                                                                           
00542   } else {                                                                      
00543     // ray interval lies completely on one side of the plane.  Test one point.
00544     float3 p = ray.origin + ray.tmin * ray.direction;                         
00545     if (dot(make_float4(p, 1.0f), plane) < 0) {
00546       ray.tmin = ray.tmax = RT_DEFAULT_MAX; // cull geometry
00547     }                                                                         
00548   }                                                                             
00549 }               
00550 
00551 
00552 //
00553 // Clear the raystats buffers to zeros
00554 //
00555 #if defined(ORT_RAYSTATS)
00556 RT_PROGRAM void clear_raystats_buffers() {
00557   raystats1_buffer[launch_index] = make_uint4(0, 0, 0, 0); // clear ray counters to zero
00558   raystats2_buffer[launch_index] = make_uint4(0, 0, 0, 0); // clear ray counters to zero
00559 }
00560 #endif
00561 
00562 
00563 //
00564 // Clear the accumulation buffer to zeros
00565 //
00566 RT_PROGRAM void clear_accumulation_buffer() {
00567   accumulation_buffer[launch_index] = make_float4(0.0f, 0.0f, 0.0f, 0.0f);
00568 }
00569 
00570 
00571 //
00572 // Copy the contents of the accumulation buffer to the destination 
00573 // framebuffer, while converting the representation of the data from
00574 // floating point down to unsigned chars, performing clamping and any
00575 // other postprocessing at the same time.
00576 //
00577 RT_PROGRAM void draw_accumulation_buffer() {
00578 #if defined(ORT_TIME_COLORING) && defined(ORT_TIME_COMBINE_MAX)
00579   float4 curcol = accumulation_buffer[launch_index];
00580 
00581   // divide time value by normalization factor (to scale up the max value)
00582   curcol.x /= accumulation_normalization_factor;
00583 
00584   // multiply the remaining color components (to average them)
00585   curcol.y *= accumulation_normalization_factor;
00586   curcol.z *= accumulation_normalization_factor;
00587   curcol.w *= accumulation_normalization_factor;
00588 
00589   framebuffer[launch_index] = make_color_rgb4u(curcol);
00590 #else
00591   framebuffer[launch_index] = make_color_rgb4u(accumulation_buffer[launch_index] * accumulation_normalization_factor);
00592 #endif
00593 }
00594 
00595 // no-op placeholder used when running with progressive rendering
00596 RT_PROGRAM void draw_accumulation_buffer_stub() {
00597 }
00598 
00599 
00600 //
00601 // OptiX programs that implement the camera models and ray generation 
00602 // code for both perspective and orthographic projection modes of VMD
00603 //
00604 
00605 //
00606 // Ray gen accumulation buffer helper routines
00607 //
00608 static void __inline__ __device__ accumulate_color(float3 &col, 
00609                                                    float alpha = 1.0f) {
00610 #if defined(VMDOPTIX_PROGRESSIVEAPI)
00611   if (progressive_enabled) {
00612     col *= accumulation_normalization_factor;
00613     alpha *= accumulation_normalization_factor;
00614 
00615 #if OPTIX_VERSION < 3080
00616     // XXX prior to OptiX 3.8, a hard-coded gamma correction was required
00617     // VCA gamma correction workaround, changes gamma 2.2 back to gamma 1.0
00618     float invgamma = 1.0f / 0.4545f;
00619     col.x = powf(col.x, invgamma);
00620     col.y = powf(col.y, invgamma);
00621     col.z = powf(col.z, invgamma);
00622 #endif
00623 
00624     // for optix-vca progressive mode accumulation is handled in server code
00625     accumulation_buffer[launch_index]  = make_float4(col, alpha);
00626   } else {
00627     // For batch mode we accumulate ourselves
00628     accumulation_buffer[launch_index] += make_float4(col, alpha);
00629   }
00630 #else
00631   // For batch mode we accumulate ourselves
00632   accumulation_buffer[launch_index] += make_float4(col, alpha);
00633 #endif
00634 }
00635 
00636 
00637 #if defined(ORT_TIME_COLORING)
00638 // special accumulation helper routine for time-based coloring
00639 static void __inline__ __device__ 
00640 accumulate_time_coloring(float3 &col, clock_t t0) {
00641   clock_t t1 = clock(); // stop per-pixel RT timer (in usec)
00642   float4 curcol = accumulation_buffer[launch_index];
00643 
00644   // overwrite the red channel with fraction of the max allowable runtime,
00645   // clamped to the range 0-1, in the case of excessively long traces
00646   float pixel_time = (t1 - t0) * ORT_TIME_NORMALIZATION;
00647   col.x = __saturatef(pixel_time);
00648 
00649 #if defined(ORT_TIME_COMBINE_MAX)
00650   // return the slowest (max) time, but average the colors
00651   curcol.x = fmaxf(curcol.x, col.x);
00652   curcol.y += col.y;
00653   curcol.z += col.z;
00654   curcol.w += 1.0f;
00655 #else
00656   // average both time and colors
00657   curcol += make_float4(col, 1.0f);
00658 #endif
00659   
00660   accumulation_buffer[launch_index] = curcol;
00661 }
00662 #endif
00663 
00664 
00665 static int __inline__ __device__ subframe_count() {
00666 #if defined(VMDOPTIX_PROGRESSIVEAPI)
00667   return (accumCount + progressiveSubframeIndex);
00668 #else
00669   return accumCount;
00670 #endif
00671 }
00672 
00673 
00674 //
00675 // CUDA device function for computing the new ray origin
00676 // and ray direction, given the radius of the circle of confusion disc,
00677 // and an orthonormal basis for each ray.
00678 //
00679 static __device__ __inline__
00680 void dof_ray(const float3 &ray_origin_orig, float3 &ray_origin, 
00681              const float3 &ray_direction_orig, float3 &ray_direction,
00682              unsigned int &randseed, const float3 &up, const float3 &right) {
00683   float3 focuspoint = ray_origin_orig + ray_direction_orig * cam_dof_focal_dist;
00684   float2 dofjxy;
00685   jitter_disc2f(randseed, dofjxy, cam_dof_aperture_rad);
00686   ray_origin = ray_origin_orig + dofjxy.x*right + dofjxy.y*up;
00687   ray_direction = normalize(focuspoint - ray_origin);
00688 }
00689 
00690 
00691 //
00692 // 360-degree stereoscopic cubemap image format for use with
00693 // Oculus, Google Cardboard, and similar VR headsets
00694 //
00695 template<int STEREO_ON, int DOF_ON> 
00696 static __device__ __inline__ 
00697 void vmd_camera_cubemap_general() {
00698 #if defined(ORT_TIME_COLORING)
00699   clock_t t0 = clock(); // start per-pixel RT timer
00700 #endif
00701 
00702   // compute which cubemap face we're drawing by the X index.
00703   uint facesz = launch_dim.y; // square cube faces, equal to image height
00704   uint face = (launch_index.x / facesz) % 6;
00705   uint2 face_idx = make_uint2(launch_index.x % facesz, launch_index.y);
00706 
00707   // For the OTOY ORBX viewer, Oculus VR software, and some of the 
00708   // related apps, the cubemap image is stored with the X axis oriented
00709   // such that when viewed as a 2-D image, they are all mirror images.
00710   // The mirrored left-right orientation used here corresponds to what is
00711   // seen standing outside the cube, whereas the ray tracer shoots
00712   // rays from the inside, so we flip the X-axis pixel storage order.
00713   // The top face of the cubemap has both the left-right and top-bottom
00714   // orientation flipped also.
00715   // Set per-face orthonormal basis for camera
00716   float3 face_U, face_V, face_W;
00717   switch (face) {
00718     case 0: // back face
00719       face_U =  cam_U;
00720       face_V =  cam_V;
00721       face_W = -cam_W;
00722       break;
00723 
00724     case 1: // front face
00725       face_U =  -cam_U;
00726       face_V =  cam_V;
00727       face_W =  cam_W;
00728       break;
00729 
00730     case 2: // top face
00731       face_U = -cam_W;
00732       face_V =  cam_U;
00733       face_W =  cam_V;
00734       break;
00735 
00736     case 3: // bottom face
00737       face_U = -cam_W;
00738       face_V = -cam_U;
00739       face_W = -cam_V;
00740       break;
00741 
00742     case 4: // left face
00743       face_U = -cam_W;
00744       face_V =  cam_V;
00745       face_W = -cam_U;
00746       break;
00747 
00748     case 5: // right face
00749       face_U =  cam_W;
00750       face_V =  cam_V;
00751       face_W =  cam_U;
00752       break;
00753   }
00754 
00755   // Stereoscopic rendering is provided by rendering in a side-by-side
00756   // format with the left eye image into the left half of a double-wide
00757   // framebuffer, and the right eye into the right half.  The subsequent 
00758   // OpenGL drawing code can trivially unpack and draw the two images 
00759   // into an efficient cubemap texture.
00760   uint viewport_sz_x, viewport_idx_x;
00761   float eyeshift;
00762   if (STEREO_ON) {
00763     // render into a double-wide framebuffer when stereo is enabled
00764     viewport_sz_x = launch_dim.x >> 1;
00765     if (launch_index.x >= viewport_sz_x) {
00766       // right image
00767       viewport_idx_x = launch_index.x - viewport_sz_x;
00768       eyeshift =  0.5f * cam_stereo_eyesep;
00769     } else {
00770       // left image
00771       viewport_idx_x = launch_index.x;
00772       eyeshift = -0.5f * cam_stereo_eyesep;
00773     }
00774   } else {
00775     // render into a normal size framebuffer if stereo is not enabled
00776     viewport_sz_x = launch_dim.x;
00777     viewport_idx_x = launch_index.x;
00778     eyeshift = 0.0f;
00779   }
00780 
00781   // 
00782   // general primary ray calculations, locked to 90-degree FoV per face...
00783   //
00784   float facescale = 1.0f / facesz;
00785   float2 d = make_float2(face_idx.x, face_idx.y) * facescale * 2.f - 1.0f; // center of pixel in image plane
00786 
00787   unsigned int randseed = tea<4>(launch_dim.x*(launch_index.y)+viewport_idx_x, subframe_count());
00788 
00789   float3 col = make_float3(0.0f);
00790   for (int s=0; s<aa_samples; s++) {
00791     float2 jxy;  
00792     jitter_offset2f(randseed, jxy);
00793     jxy = jxy * facescale * 2.f + d;
00794     float3 ray_direction = normalize(jxy.x*face_U + jxy.y*face_V + face_W);
00795 
00796     float3 ray_origin = cam_pos;
00797     if (STEREO_ON) {
00798       ray_origin += eyeshift * cross(ray_direction, cam_V);
00799     }
00800 
00801     // compute new ray origin and ray direction
00802     if (DOF_ON) {
00803       dof_ray(ray_origin, ray_origin, ray_direction, ray_direction,
00804               randseed, face_V, face_U);
00805     }
00806 
00807     // trace the new ray...
00808     PerRayData_radiance prd;
00809     prd.importance = 1.f;
00810     prd.depth = 0;
00811     prd.transcnt = max_trans;
00812     optix::Ray ray = optix::make_Ray(ray_origin, ray_direction, radiance_ray_type, scene_epsilon, RT_DEFAULT_MAX);
00813     rtTrace(root_object, ray, prd);
00814     col += prd.result; 
00815   }
00816 #if defined(ORT_RAYSTATS)
00817   raystats1_buffer[launch_index].x+=aa_samples; // increment primary ray counter
00818 #endif
00819 
00820 #if defined(ORT_TIME_COLORING)
00821   accumulate_time_coloring(col, t0);
00822 #else
00823   accumulate_color(col);
00824 #endif
00825 }
00826 
00827 RT_PROGRAM void vmd_camera_cubemap() {
00828   vmd_camera_cubemap_general<0, 0>();
00829 }
00830 
00831 RT_PROGRAM void vmd_camera_cubemap_dof() {
00832   vmd_camera_cubemap_general<0, 1>();
00833 }
00834 
00835 RT_PROGRAM void vmd_camera_cubemap_stereo() {
00836   vmd_camera_cubemap_general<1, 0>();
00837 }
00838 
00839 RT_PROGRAM void vmd_camera_cubemap_stereo_dof() {
00840   vmd_camera_cubemap_general<1, 1>();
00841 }
00842 
00843 
00844 
00845 //
00846 // Camera ray generation code for planetarium dome display
00847 // Generates a fisheye style frame with ~180 degree FoV
00848 //
00849 template<int STEREO_ON, int DOF_ON>
00850 static __device__ __inline__
00851 void vmd_camera_dome_general() {
00852 #if defined(ORT_TIME_COLORING)
00853   clock_t t0 = clock(); // start per-pixel RT timer
00854 #endif
00855 
00856   // Stereoscopic rendering is provided by rendering in an over/under
00857   // format with the left eye image into the top half of a double-high
00858   // framebuffer, and the right eye into the lower half.  The subsequent
00859   // OpenGL drawing code can trivially unpack and draw the two images
00860   // with simple pointer offset arithmetic.
00861   uint viewport_sz_y, viewport_idx_y;
00862   float eyeshift;
00863   if (STEREO_ON) {
00864     // render into a double-high framebuffer when stereo is enabled
00865     viewport_sz_y = launch_dim.y >> 1;
00866     if (launch_index.y >= viewport_sz_y) {
00867       // left image
00868       viewport_idx_y = launch_index.y - viewport_sz_y;
00869       eyeshift = -0.5f * cam_stereo_eyesep;
00870     } else {
00871       // right image
00872       viewport_idx_y = launch_index.y;
00873       eyeshift =  0.5f * cam_stereo_eyesep;
00874     }
00875   } else {
00876     // render into a normal size framebuffer if stereo is not enabled
00877     viewport_sz_y = launch_dim.y;
00878     viewport_idx_y = launch_index.y;
00879     eyeshift = 0.0f;
00880   }
00881 
00882   float fov = M_PIf; // dome FoV in radians
00883 
00884   // half FoV in radians, pixels beyond this distance are outside
00885   // of the field of view of the projection, and are set black
00886   float thetamax = 0.5 * fov;
00887 
00888   // The dome angle from center of the projection is proportional
00889   // to the image-space distance from the center of the viewport.
00890   // viewport_sz contains the viewport size, radperpix contains the
00891   // radians/pixel scaling factors in X/Y, and viewport_mid contains
00892   // the midpoint coordinate of the viewpoint used to compute the
00893   // distance from center.
00894   float2 viewport_sz = make_float2(launch_dim.x, viewport_sz_y);
00895   float2 radperpix = fov / viewport_sz;
00896   float2 viewport_mid = viewport_sz * 0.5f;
00897 
00898   unsigned int randseed = tea<4>(launch_dim.x*(launch_index.y)+launch_index.x, subframe_count());
00899 
00900   float3 col = make_float3(0.0f);
00901   float alpha = 0.0f;
00902   for (int s=0; s<aa_samples; s++) {
00903     // compute the jittered image plane sample coordinate
00904     float2 jxy;
00905     jitter_offset2f(randseed, jxy);
00906     float2 viewport_idx = make_float2(launch_index.x, viewport_idx_y) + jxy;
00907 
00908     // compute the ray angles in X/Y and total angular distance from center
00909     float2 p = (viewport_idx - viewport_mid) * radperpix;
00910     float theta = hypotf(p.x, p.y);
00911 
00912     // pixels outside the dome FoV are treated as black by not
00913     // contributing to the color accumulator
00914     if (theta < thetamax) {
00915       float3 ray_direction;
00916       float3 ray_origin = cam_pos;
00917 
00918       if (theta == 0) {
00919         // handle center of dome where azimuth is undefined by
00920         // setting the ray direction to the zenith
00921         ray_direction = cam_W;
00922       } else {
00923         float sintheta, costheta;
00924         sincosf(theta, &sintheta, &costheta);
00925         float rsin = sintheta / theta; // normalize component
00926         ray_direction = cam_U*rsin*p.x + cam_V*rsin*p.y + cam_W*costheta;
00927         if (STEREO_ON) {
00928           // assumes a flat dome, where cam_W also points in the
00929           // audience "up" direction
00930           ray_origin += eyeshift * cross(ray_direction, cam_W);
00931         }
00932 
00933         if (DOF_ON) {
00934           float rcos = costheta / theta; // normalize component
00935           float3 ray_up    = -cam_U*rcos*p.x  -cam_V*rcos*p.y + cam_W*sintheta;
00936           float3 ray_right =  cam_U*(p.y/theta) + cam_V*(-p.x/theta);
00937           dof_ray(ray_origin, ray_origin, ray_direction, ray_direction,
00938                   randseed, ray_up, ray_right);
00939         }
00940       }
00941 
00942       // trace the new ray...
00943       PerRayData_radiance prd;
00944       prd.importance = 1.f;
00945       prd.alpha = 1.f;
00946       prd.depth = 0;
00947       prd.transcnt = max_trans;
00948       optix::Ray ray = optix::make_Ray(ray_origin, ray_direction, radiance_ray_type, scene_epsilon, RT_DEFAULT_MAX);
00949       rtTrace(root_object, ray, prd);
00950       col += prd.result;
00951       alpha += prd.alpha;
00952     }
00953   }
00954 #if defined(ORT_RAYSTATS)
00955   raystats1_buffer[launch_index].x+=aa_samples; // increment primary ray counter
00956 #endif
00957 
00958 #if defined(ORT_TIME_COLORING)
00959   accumulate_time_coloring(col, t0);
00960 #else
00961   accumulate_color(col, alpha);
00962 #endif
00963 }
00964 
00965 RT_PROGRAM void vmd_camera_dome_master() {
00966   vmd_camera_dome_general<0, 0>();
00967 }
00968 
00969 RT_PROGRAM void vmd_camera_dome_master_dof() {
00970   vmd_camera_dome_general<0, 1>();
00971 }
00972 
00973 RT_PROGRAM void vmd_camera_dome_master_stereo() {
00974   vmd_camera_dome_general<1, 0>();
00975 }
00976 
00977 RT_PROGRAM void vmd_camera_dome_master_stereo_dof() {
00978   vmd_camera_dome_general<1, 1>();
00979 }
00980 
00981 
00982 //
00983 // Camera ray generation code for 360 degre FoV 
00984 // equirectangular (lat/long) projection suitable
00985 // for use a texture map for a sphere, e.g. for 
00986 // immersive VR HMDs, other spheremap-based projections.
00987 // 
00988 template<int STEREO_ON, int DOF_ON>
00989 static __device__ __inline__
00990 void vmd_camera_equirectangular_general() {
00991 #if defined(ORT_TIME_COLORING)
00992   clock_t t0 = clock(); // start per-pixel RT timer
00993 #endif
00994 
00995   // The Samsung GearVR OTOY ORBX players have the left eye image on top, 
00996   // and the right eye image on the bottom.
00997   // Stereoscopic rendering is provided by rendering in an over/under
00998   // format with the left eye image into the top half of a double-high 
00999   // framebuffer, and the right eye into the lower half.  The subsequent 
01000   // OpenGL drawing code can trivially unpack and draw the two images 
01001   // with simple pointer offset arithmetic.
01002   uint viewport_sz_y, viewport_idx_y;
01003   float eyeshift;
01004   if (STEREO_ON) {
01005     // render into a double-high framebuffer when stereo is enabled
01006     viewport_sz_y = launch_dim.y >> 1;
01007     if (launch_index.y >= viewport_sz_y) {
01008       // left image
01009       viewport_idx_y = launch_index.y - viewport_sz_y;
01010       eyeshift = -0.5f * cam_stereo_eyesep;
01011     } else {
01012       // right image
01013       viewport_idx_y = launch_index.y;
01014       eyeshift =  0.5f * cam_stereo_eyesep;
01015     }
01016   } else {
01017     // render into a normal size framebuffer if stereo is not enabled
01018     viewport_sz_y = launch_dim.y;
01019     viewport_idx_y = launch_index.y;
01020     eyeshift = 0.0f;
01021   }
01022 
01023   float2 viewport_sz = make_float2(launch_dim.x, viewport_sz_y);
01024   float2 radperpix = M_PIf / viewport_sz * make_float2(2.0f, 1.0f);
01025   float2 viewport_mid = viewport_sz * 0.5f;
01026 
01027   unsigned int randseed = tea<4>(launch_dim.x*(launch_index.y)+launch_index.x, subframe_count());
01028 
01029   float3 col = make_float3(0.0f);
01030   for (int s=0; s<aa_samples; s++) {
01031     float2 jxy;  
01032     jitter_offset2f(randseed, jxy);
01033 
01034     float2 viewport_idx = make_float2(launch_index.x, viewport_idx_y) + jxy;
01035     float2 rangle = (viewport_idx - viewport_mid) * radperpix;
01036 
01037     float sin_ax, cos_ax, sin_ay, cos_ay;
01038     sincosf(rangle.x, &sin_ax, &cos_ax);
01039     sincosf(rangle.y, &sin_ay, &cos_ay);
01040 
01041     float3 ray_direction = normalize(cos_ay * (cos_ax * cam_W + sin_ax * cam_U) + sin_ay * cam_V);
01042 
01043     float3 ray_origin = cam_pos;
01044     if (STEREO_ON) {
01045       ray_origin += eyeshift * cross(ray_direction, cam_V);
01046     }
01047 
01048     // compute new ray origin and ray direction
01049     if (DOF_ON) {
01050       float3 ray_right = normalize(cos_ay * (-sin_ax * cam_W - cos_ax * cam_U) + sin_ay * cam_V);
01051       float3 ray_up = cross(ray_direction, ray_right);
01052       dof_ray(ray_origin, ray_origin, ray_direction, ray_direction,
01053               randseed, ray_up, ray_right);
01054     }
01055 
01056     // trace the new ray...
01057     PerRayData_radiance prd;
01058     prd.importance = 1.f;
01059     prd.depth = 0;
01060     prd.transcnt = max_trans;
01061     optix::Ray ray = optix::make_Ray(ray_origin, ray_direction, radiance_ray_type, scene_epsilon, RT_DEFAULT_MAX);
01062     rtTrace(root_object, ray, prd);
01063     col += prd.result;
01064   }
01065 #if defined(ORT_RAYSTATS)
01066   raystats1_buffer[launch_index].x+=aa_samples; // increment primary ray counter
01067 #endif
01068 
01069 #if defined(ORT_TIME_COLORING)
01070   accumulate_time_coloring(col, t0);
01071 #else
01072   accumulate_color(col);
01073 #endif
01074 }
01075 
01076 RT_PROGRAM void vmd_camera_equirectangular() {
01077   vmd_camera_equirectangular_general<0, 0>();
01078 }
01079 
01080 RT_PROGRAM void vmd_camera_equirectangular_dof() {
01081   vmd_camera_equirectangular_general<0, 1>();
01082 }
01083 
01084 RT_PROGRAM void vmd_camera_equirectangular_stereo() {
01085   vmd_camera_equirectangular_general<1, 0>();
01086 }
01087 
01088 RT_PROGRAM void vmd_camera_equirectangular_stereo_dof() {
01089   vmd_camera_equirectangular_general<1, 1>();
01090 }
01091 
01092 
01093 //
01094 // Templated Oculus Rift perspective camera ray generation code
01095 //
01096 template<int STEREO_ON, int DOF_ON> 
01097 static __device__ __inline__ 
01098 void vmd_camera_oculus_rift_general() {
01099 #if defined(ORT_TIME_COLORING)
01100   clock_t t0 = clock(); // start per-pixel RT timer
01101 #endif
01102 
01103   // Stereoscopic rendering is provided by rendering in a side-by-side
01104   // format with the left eye image in the left half of a double-wide
01105   // framebuffer, and the right eye in the right half.  The subsequent 
01106   // OpenGL drawing code can trivially unpack and draw the two images 
01107   // with simple pointer offset arithmetic.
01108   uint viewport_sz_x, viewport_idx_x;
01109   float eyeshift;
01110   if (STEREO_ON) {
01111     // render into a double-wide framebuffer when stereo is enabled
01112     viewport_sz_x = launch_dim.x >> 1;
01113     if (launch_index.x >= viewport_sz_x) {
01114       // right image
01115       viewport_idx_x = launch_index.x - viewport_sz_x;
01116       eyeshift =  0.5f * cam_stereo_eyesep;
01117     } else {
01118       // left image
01119       viewport_idx_x = launch_index.x;
01120       eyeshift = -0.5f * cam_stereo_eyesep;
01121     }
01122   } else {
01123     // render into a normal size framebuffer if stereo is not enabled
01124     viewport_sz_x = launch_dim.x;
01125     viewport_idx_x = launch_index.x;
01126     eyeshift = 0.0f;
01127   }
01128 
01129   // 
01130   // general primary ray calculations
01131   //
01132   float2 aspect = make_float2(float(viewport_sz_x) / float(launch_dim.y), 1.0f) * cam_zoom;
01133   float2 viewportscale = 1.0f / make_float2(viewport_sz_x, launch_dim.y);
01134   float2 d = make_float2(viewport_idx_x, launch_index.y) * viewportscale * aspect * 2.f - aspect; // center of pixel in image plane
01135 
01136 
01137   // Compute barrel distortion required to correct for the pincushion inherent
01138   // in the plano-convex optics in the Oculus Rift, Google Cardboard, etc.
01139   // Barrel distortion involves computing distance of the pixel from the 
01140   // center of the eye viewport, and then scaling this distance by a factor
01141   // based on the original distance: 
01142   //   rnew = 0.24 * r^4 + 0.22 * r^2 + 1.0
01143   // Since we are only using even powers of r, we can use efficient 
01144   // squared distances everywhere.
01145   // The current implementation doesn't discard rays that would have fallen
01146   // outside of the original viewport FoV like most OpenGL implementations do.
01147   // The current implementation computes the distortion for the initial ray 
01148   // but doesn't apply these same corrections to antialiasing jitter, to
01149   // depth-of-field jitter, etc, so this leaves something to be desired if
01150   // we want best quality, but this raygen code is really intended for 
01151   // interactive display on an Oculus Rift or Google Cardboard type viewer,
01152   // so I err on the side of simplicity/speed for now. 
01153   float2 cp = make_float2(viewport_sz_x >> 1, launch_dim.y >> 1) * viewportscale * aspect * 2.f - aspect;;
01154   float2 dr = d - cp;
01155   float r2 = dr.x*dr.x + dr.y*dr.y;
01156   float r = 0.24f*r2*r2 + 0.22f*r2 + 1.0f;
01157   d = r * dr; 
01158 
01159   int subframecount = subframe_count();
01160   unsigned int randseed = tea<4>(launch_dim.x*(launch_index.y)+viewport_idx_x, subframecount);
01161 
01162   float3 eyepos = cam_pos;
01163   if (STEREO_ON) {
01164     eyepos += eyeshift * cam_U;
01165   } 
01166 
01167   float3 ray_origin = eyepos;
01168   float3 col = make_float3(0.0f);
01169   for (int s=0; s<aa_samples; s++) {
01170     float2 jxy;  
01171     jitter_offset2f(randseed, jxy);
01172 
01173     // don't jitter the first sample, since when using an HMD we often run
01174     // with only one sample per pixel unless the user wants higher fidelity
01175     jxy *= (subframecount > 0 || s > 0);
01176 
01177     jxy = jxy * viewportscale * aspect * 2.f + d;
01178     float3 ray_direction = normalize(jxy.x*cam_U + jxy.y*cam_V + cam_W);
01179  
01180     // compute new ray origin and ray direction
01181     if (DOF_ON) {
01182       dof_ray(eyepos, ray_origin, ray_direction, ray_direction,
01183               randseed, cam_V, cam_U);
01184     }
01185 
01186     // trace the new ray...
01187     PerRayData_radiance prd;
01188     prd.importance = 1.f;
01189     prd.depth = 0;
01190     prd.transcnt = max_trans;
01191     optix::Ray ray = optix::make_Ray(ray_origin, ray_direction, radiance_ray_type, scene_epsilon, RT_DEFAULT_MAX);
01192     rtTrace(root_object, ray, prd);
01193     col += prd.result; 
01194   }
01195 #if defined(ORT_RAYSTATS)
01196   raystats1_buffer[launch_index].x+=aa_samples; // increment primary ray counter
01197 #endif
01198 
01199 #if defined(ORT_TIME_COLORING)
01200   accumulate_time_coloring(col, t0);
01201 #else
01202   accumulate_color(col);
01203 #endif
01204 }
01205 
01206 RT_PROGRAM void vmd_camera_oculus_rift() {
01207   vmd_camera_oculus_rift_general<0, 0>();
01208 }
01209 
01210 RT_PROGRAM void vmd_camera_oculus_rift_dof() {
01211   vmd_camera_oculus_rift_general<0, 1>();
01212 }
01213 
01214 RT_PROGRAM void vmd_camera_oculus_rift_stereo() {
01215   vmd_camera_oculus_rift_general<1, 0>();
01216 }
01217 
01218 RT_PROGRAM void vmd_camera_oculus_rift_stereo_dof() {
01219   vmd_camera_oculus_rift_general<1, 1>();
01220 }
01221 
01222 
01223 
01224 //
01225 // Templated perspective camera ray generation code
01226 //
01227 template<int STEREO_ON, int DOF_ON> 
01228 static __device__ __inline__ 
01229 void vmd_camera_perspective_general() {
01230 #if defined(ORT_TIME_COLORING)
01231   clock_t t0 = clock(); // start per-pixel RT timer
01232 #endif
01233 
01234   // Stereoscopic rendering is provided by rendering in an over/under
01235   // format with the left eye image into the top half of a double-high 
01236   // framebuffer, and the right eye into the lower half.  The subsequent 
01237   // OpenGL drawing code can trivially unpack and draw the two images 
01238   // with simple pointer offset arithmetic.
01239   float3 eyepos;
01240   uint viewport_sz_y, viewport_idx_y;
01241   if (STEREO_ON) {
01242     // render into a double-high framebuffer when stereo is enabled
01243     viewport_sz_y = launch_dim.y >> 1;
01244     if (launch_index.y >= viewport_sz_y) {
01245       // right image
01246       viewport_idx_y = launch_index.y - viewport_sz_y;
01247       eyepos = cam_pos + cam_U * cam_stereo_eyesep * 0.5f;
01248     } else {
01249       // left image
01250       viewport_idx_y = launch_index.y;
01251       eyepos = cam_pos - cam_U * cam_stereo_eyesep * 0.5f;
01252     }
01253   } else {
01254     // render into a normal size framebuffer if stereo is not enabled
01255     viewport_sz_y = launch_dim.y;
01256     viewport_idx_y = launch_index.y;
01257     eyepos = cam_pos;
01258   }
01259 
01260   // 
01261   // general primary ray calculations
01262   //
01263   float2 aspect = make_float2(float(launch_dim.x) / float(viewport_sz_y), 1.0f) * cam_zoom;
01264   float2 viewportscale = 1.0f / make_float2(launch_dim.x, viewport_sz_y);
01265   float2 d = make_float2(launch_index.x, viewport_idx_y) * viewportscale * aspect * 2.f - aspect; // center of pixel in image plane
01266 
01267   unsigned int randseed = tea<4>(launch_dim.x*(viewport_idx_y)+launch_index.x, subframe_count());
01268 
01269   float3 col = make_float3(0.0f);
01270   float alpha = 0.0f;
01271   float3 ray_origin = eyepos;
01272   for (int s=0; s<aa_samples; s++) {
01273     float2 jxy;  
01274     jitter_offset2f(randseed, jxy);
01275 
01276     jxy = jxy * viewportscale * aspect * 2.f + d;
01277     float3 ray_direction = normalize(jxy.x*cam_U + jxy.y*cam_V + cam_W);
01278 
01279     // compute new ray origin and ray direction
01280     if (DOF_ON) {
01281       dof_ray(eyepos, ray_origin, ray_direction, ray_direction,
01282               randseed, cam_V, cam_U);
01283     }
01284 
01285     // trace the new ray...
01286     PerRayData_radiance prd;
01287     prd.importance = 1.f;
01288     prd.alpha = 1.f;
01289     prd.depth = 0;
01290     prd.transcnt = max_trans;
01291     optix::Ray ray = optix::make_Ray(ray_origin, ray_direction, radiance_ray_type, scene_epsilon, RT_DEFAULT_MAX);
01292     rtTrace(root_object, ray, prd);
01293     col += prd.result; 
01294     alpha += prd.alpha;
01295   }
01296 #if defined(ORT_RAYSTATS)
01297   raystats1_buffer[launch_index].x+=aa_samples; // increment primary ray counter
01298 #endif
01299 
01300 #if defined(ORT_TIME_COLORING)
01301   accumulate_time_coloring(col, t0);
01302 #else
01303   accumulate_color(col, alpha);
01304 #endif
01305 }
01306 
01307 
01308 RT_PROGRAM void vmd_camera_perspective() {
01309   vmd_camera_perspective_general<0, 0>();
01310 }
01311 
01312 RT_PROGRAM void vmd_camera_perspective_dof() {
01313   vmd_camera_perspective_general<0, 1>();
01314 }
01315 
01316 RT_PROGRAM void vmd_camera_perspective_stereo() {
01317   vmd_camera_perspective_general<1, 0>();
01318 }
01319 
01320 RT_PROGRAM void vmd_camera_perspective_stereo_dof() {
01321   vmd_camera_perspective_general<1, 1>();
01322 }
01323 
01324 
01325 //
01326 // Templated orthographic camera ray generation code
01327 //
01328 template<int STEREO_ON, int DOF_ON> 
01329 static __device__ __inline__ 
01330 void vmd_camera_orthographic_general() {
01331 #if defined(ORT_TIME_COLORING)
01332   clock_t t0 = clock(); // start per-pixel RT timer
01333 #endif
01334 
01335   // Stereoscopic rendering is provided by rendering in an over/under
01336   // format with the left eye image into the top half of a double-high 
01337   // framebuffer, and the right eye into the lower half.  The subsequent 
01338   // OpenGL drawing code can trivially unpack and draw the two images 
01339   // with simple pointer offset arithmetic.
01340   float3 eyepos;
01341   uint viewport_sz_y, viewport_idx_y;
01342   float3 view_direction;
01343   if (STEREO_ON) {
01344     // render into a double-high framebuffer when stereo is enabled
01345     viewport_sz_y = launch_dim.y >> 1;
01346     if (launch_index.y >= viewport_sz_y) {
01347       // right image
01348       viewport_idx_y = launch_index.y - viewport_sz_y;
01349       eyepos = cam_pos + cam_U * cam_stereo_eyesep * 0.5f;
01350     } else {
01351       // left image
01352       viewport_idx_y = launch_index.y;
01353       eyepos = cam_pos - cam_U * cam_stereo_eyesep * 0.5f;
01354     }
01355     view_direction = normalize(cam_pos-eyepos + normalize(cam_W) * cam_stereo_convergence_dist);
01356   } else {
01357     // render into a normal size framebuffer if stereo is not enabled
01358     viewport_sz_y = launch_dim.y;
01359     viewport_idx_y = launch_index.y;
01360     eyepos = cam_pos;
01361     view_direction = normalize(cam_W);
01362   }
01363 
01364   // 
01365   // general primary ray calculations
01366   //
01367   float2 aspect = make_float2(float(launch_dim.x) / float(viewport_sz_y), 1.0f) * cam_zoom;
01368   float2 viewportscale = 1.0f / make_float2(launch_dim.x, viewport_sz_y);
01369 
01370   float2 d = make_float2(launch_index.x, viewport_idx_y) * viewportscale * aspect * 2.f - aspect; // center of pixel in image plane
01371 
01372   unsigned int randseed = tea<4>(launch_dim.x*(viewport_idx_y)+launch_index.x, subframe_count());
01373 
01374   float3 col = make_float3(0.0f);
01375   float alpha = 0.0f;
01376   float3 ray_direction = view_direction;
01377   for (int s=0; s<aa_samples; s++) {
01378     float2 jxy;  
01379     jitter_offset2f(randseed, jxy);
01380     jxy = jxy * viewportscale * aspect * 2.f + d;
01381     float3 ray_origin = eyepos + jxy.x*cam_U + jxy.y*cam_V;
01382 
01383     // compute new ray origin and ray direction
01384     if (DOF_ON) {
01385       dof_ray(ray_origin, ray_origin, view_direction, ray_direction,
01386               randseed, cam_V, cam_U);
01387     }
01388 
01389     // trace the new ray...
01390     PerRayData_radiance prd;
01391     prd.importance = 1.f;
01392     prd.alpha = 1.f;
01393     prd.depth = 0;
01394     prd.transcnt = max_trans;
01395     optix::Ray ray = optix::make_Ray(ray_origin, ray_direction, radiance_ray_type, scene_epsilon, RT_DEFAULT_MAX);
01396     rtTrace(root_object, ray, prd);
01397     col += prd.result; 
01398     alpha += prd.alpha;
01399   }
01400 #if defined(ORT_RAYSTATS)
01401   raystats1_buffer[launch_index].x+=aa_samples; // increment primary ray counter
01402 #endif
01403 
01404 #if defined(ORT_TIME_COLORING)
01405   accumulate_time_coloring(col, t0);
01406 #else
01407   accumulate_color(col, alpha);
01408 #endif
01409 }
01410 
01411 RT_PROGRAM void vmd_camera_orthographic() {
01412   vmd_camera_orthographic_general<0, 0>();
01413 }
01414 
01415 RT_PROGRAM void vmd_camera_orthographic_dof() {
01416   vmd_camera_orthographic_general<0, 1>();
01417 }
01418 
01419 RT_PROGRAM void vmd_camera_orthographic_stereo() {
01420   vmd_camera_orthographic_general<1, 0>();
01421 }
01422 
01423 RT_PROGRAM void vmd_camera_orthographic_stereo_dof() {
01424   vmd_camera_orthographic_general<1, 1>();
01425 }
01426 
01427 
01428 //
01429 // Default exception handling behavior
01430 //
01431 RT_PROGRAM void exception() {
01432   const unsigned int code = rtGetExceptionCode();
01433   switch (code) {
01434     case RT_EXCEPTION_STACK_OVERFLOW:
01435       rtPrintf("Stack overflow at launch index (%d,%d):\n",
01436                launch_index.x, launch_index.y );
01437       break;
01438 
01439 #if OPTIX_VERSION >= 3050
01440     case RT_EXCEPTION_TEXTURE_ID_INVALID:
01441     case RT_EXCEPTION_BUFFER_ID_INVALID:
01442 #endif
01443     case RT_EXCEPTION_INDEX_OUT_OF_BOUNDS:
01444     case RT_EXCEPTION_BUFFER_INDEX_OUT_OF_BOUNDS:
01445     case RT_EXCEPTION_INVALID_RAY:
01446     case RT_EXCEPTION_INTERNAL_ERROR:
01447     case RT_EXCEPTION_USER:
01448     default: 
01449       printf("Caught exception 0x%X (%d) at launch index (%d,%d)\n", 
01450              code, code, launch_index.x, launch_index.y );
01451       break;
01452   }
01453 
01454 #ifndef VMDOPTIX_PROGRESSIVEAPI
01455   framebuffer[launch_index] = make_color_rgb4u(make_float3(0.f, 0.f, 0.f));
01456 #endif
01457 }
01458 
01459 
01460 //
01461 // Fog helper function based on Tachyon
01462 //
01463 static __device__ float fog_coord(float3 hit_point) {
01464 #if defined(CADENSMOVIE)
01465   // For planetarium dome renderings, use spherical distance rather than planar
01466   // XXX this should be a runtime parameter eventually
01467   float r = t_hit;
01468 #else
01469   // Compute planar fog (e.g. to match OpenGL) by projecting t value onto 
01470   // the camera view direction vector to yield a planar a depth value.
01471   float r = dot(ray.direction, cam_W) * t_hit;
01472 #endif
01473 
01474   float f=1.0f;
01475   float v;
01476 
01477   switch (fog_mode) { 
01478     case 1: // RT_FOG_LINEAR
01479       f = (fog_end - r) / (fog_end - fog_start);
01480       break;
01481 
01482     case 2: // RT_FOG_EXP
01483       // XXX Tachyon allows fog_start to be non-zero for exponential fog,
01484       //     but OpenGL and VMD do not...
01485       // float v = fog_density * (r - fog_start);
01486       v = fog_density * r;
01487       f = expf(-v);
01488       break;
01489 
01490     case 3: // RT_FOG_EXP2
01491       // XXX Tachyon allows fog_start to be non-zero for exponential fog,
01492       //     but OpenGL and VMD do not...
01493       // float v = fog_density * (r - fog_start);
01494       v = fog_density * r;
01495       f = expf(-v*v);
01496       break;
01497 
01498     case 0: // RT_FOG_NONE
01499     default:
01500       break;
01501   }
01502 
01503   return __saturatef(f);
01504 }
01505 
01506 
01507 static __device__ float3 fog_color(float fogmod, float3 hit_col) {
01508   float3 col = (fogmod * hit_col) + ((1.0f - fogmod) * scene_bg_color);
01509   return col;
01510 }
01511 
01512 
01513 
01514 //
01515 // trivial ambient occlusion implementation based on Tachyon
01516 //
01517 static __device__ float3 shade_ambient_occlusion(float3 hit, float3 N, float aoimportance) {
01518   // unweighted non-importance-sampled scaling factor
01519   float lightscale = 2.0f / ao_samples;
01520   float3 inten = make_float3(0.0f);
01521 
01522   unsigned int randseed = tea<2>(subframe_count(), subframe_count()); 
01523 
01524   PerRayData_shadow shadow_prd;
01525 #if 1
01526   // do all the samples requested, with no observance of importance
01527   for (int s=0; s<ao_samples; s++) {
01528 #else
01529   // dynamically scale the number of AO samples depending on the 
01530   // importance assigned to the incident ray and the opacity of the 
01531   // surface we are lighting.
01532   // XXX this scheme can create crawlies when animating since the 
01533   //     AO sample rays are no longer identical between neighboring 
01534   //     pixels and there's no guarantee that the samples we're skipping
01535   //     were low-importance in terms of their contribution.
01536   //     This kind of scheme would need much more development to be usable.
01537   int nsamples = ao_samples * prd.importance * aoimportance;
01538   if (nsamples < 1)
01539     nsamples=1;
01540   lightscale = 2.0f / nsamples;
01541   for (int s=0; s<nsamples; s++) {
01542 #endif
01543     float3 dir;
01544     jitter_sphere3f(randseed, dir);
01545     float ndotambl = dot(N, dir);
01546 
01547     // flip the ray so it's in the same hemisphere as the surface normal
01548     if (ndotambl < 0.0f) {
01549       ndotambl = -ndotambl;
01550       dir = -dir;
01551     }
01552 
01553     Ray ambray;
01554 #ifdef USE_REVERSE_SHADOW_RAYS 
01555     if (shadows_enabled == RT_SHADOWS_ON_REVERSE) {
01556       // reverse any-hit ray traversal direction for increased perf
01557       // XXX We currently hard-code REVERSE_RAY_LENGTH in such a way that
01558       //     it works well for scenes that fall within the VMD view volume,
01559       //     given the relationship between the model and camera coordinate
01560       //     systems, but this would be best computed by the diagonal of the 
01561       //     AABB for the full scene, and then scaled into camera coordinates.
01562       //     The REVERSE_RAY_STEP size is computed to avoid self intersection 
01563       //     with the surface we're shading.
01564       float tmax = REVERSE_RAY_LENGTH - REVERSE_RAY_STEP;
01565       ambray = make_Ray(hit + dir * REVERSE_RAY_LENGTH, -dir, shadow_ray_type, 0, tmax);
01566     } else
01567 #endif
01568 #if defined(ORT_USE_RAY_STEP) 
01569     ambray = make_Ray(hit + ORT_RAY_STEP, dir, shadow_ray_type, scene_epsilon, ao_maxdist);
01570 #else
01571     ambray = make_Ray(hit, dir, shadow_ray_type, scene_epsilon, ao_maxdist);
01572 #endif
01573     shadow_prd.attenuation = make_float3(1.0f);
01574 
01575     // NOTE: There are no safety checks for maximum RT recursion depth
01576     //       here, so we must have asked OptiX for a max recursion depth
01577     //       that is one greater than we're tracking elsewhere internally.
01578     rtTrace(root_shadower, ambray, shadow_prd);
01579     inten += ndotambl * shadow_prd.attenuation;
01580   } 
01581 #if defined(ORT_RAYSTATS)
01582   raystats1_buffer[launch_index].z+=ao_samples; // increment AO shadow ray counter
01583 #endif
01584 
01585   return inten * lightscale;
01586 }
01587 
01588 
01589 template<int SHADOWS_ON>       
01590 static __device__ __inline__ void shade_light(float3 &result,
01591                                               float3 &hit_point, 
01592                                               float3 &N, float3 &L, 
01593                                               float p_Kd, 
01594                                               float p_Ks,
01595                                               float p_phong_exp,
01596                                               float3 &col, 
01597                                               float3 &phongcol,
01598                                               float shadow_tmax) {
01599   float inten = dot(N, L);
01600 
01601   // cast shadow ray
01602   float3 light_attenuation = make_float3(static_cast<float>(inten > 0.0f));
01603   if (SHADOWS_ON && shadows_enabled && inten > 0.0f) {
01604     PerRayData_shadow shadow_prd;
01605     shadow_prd.attenuation = make_float3(1.0f);
01606 
01607     Ray shadow_ray;
01608 #ifdef USE_REVERSE_SHADOW_RAYS
01609     if (shadows_enabled == RT_SHADOWS_ON_REVERSE) {
01610       // reverse any-hit ray traversal direction for increased perf
01611       // XXX We currently hard-code REVERSE_RAY_LENGTH in such a way that
01612       //     it works well for scenes that fall within the VMD view volume,
01613       //     given the relationship between the model and camera coordinate
01614       //     systems, but this would be best computed by the diagonal of the 
01615       //     AABB for the full scene, and then scaled into camera coordinates.
01616       //     The REVERSE_RAY_STEP size is computed to avoid self intersection 
01617       //     with the surface we're shading.
01618       float tmax = REVERSE_RAY_LENGTH - REVERSE_RAY_STEP;
01619       shadow_ray = make_Ray(hit_point + L * REVERSE_RAY_LENGTH, -L, shadow_ray_type, 0, fminf(tmax, shadow_tmax)); 
01620     } 
01621     else
01622 #endif
01623     shadow_ray = make_Ray(hit_point + ORT_RAY_STEP, L, shadow_ray_type, scene_epsilon, shadow_tmax);
01624 
01625     // NOTE: There are no safety checks for maximum RT recursion depth
01626     //       here, so we must have asked OptiX for a max recursion depth
01627     //       that is one greater than we're tracking elsewhere internally.
01628     rtTrace(root_shadower, shadow_ray, shadow_prd);
01629 #if defined(ORT_RAYSTATS)
01630     raystats1_buffer[launch_index].y++; // increment shadow ray counter
01631 #endif
01632     light_attenuation = shadow_prd.attenuation;
01633   }
01634 
01635   // If not completely shadowed, light the hit point.
01636   // When shadows are disabled, the light can't possibly be attenuated.
01637   if (!SHADOWS_ON || fmaxf(light_attenuation) > 0.0f) {
01638     result += col * p_Kd * inten * light_attenuation;
01639 
01640     // add specular hightlight using Blinn's halfway vector approach
01641     float3 H = normalize(L - ray.direction);
01642     float nDh = dot(N, H);
01643     if (nDh > 0) {
01644       float power = powf(nDh, p_phong_exp);
01645       phongcol += make_float3(p_Ks) * power * light_attenuation;
01646     }
01647   }
01648 }
01649 
01650 
01651 //
01652 // Minimalistic VMD-centric re-implementation of the key portions 
01653 // of Tachyon's main shader.  
01654 //
01655 // This shader has been written to be expanded into a large set of
01656 // fully specialized shaders generated through combinatorial expansion
01657 // of each of the major shader features associated with scene-wide or
01658 // material-specific shading properties.
01659 // At present, there are three scene-wide properties (fog, shadows, AO),
01660 // and three material-specific properties (outline, reflection, transmission).
01661 // Tere can be a performance cost for OptiX work scheduling of disparate 
01662 // materials if too many unique materials are used in a scene. 
01663 // Although there are 8 combinations of scene-wide parameters and 
01664 // 8 combinations of material-specific parameters (64 in total),
01665 // the scene-wide parameters are uniform for the whole scene. 
01666 // We will therefore only have at most 8 different shader variants 
01667 // in use in a given scene, due to the 8 possible combinations
01668 // of material-specific (outline, reflection, transmission) properties.
01669 // 
01670 // The macros that generate the full set of 64 possible shader variants
01671 // are at the very end of this source file.
01672 //
01673 template<int CLIP_VIEW_ON,     
01674          int HEADLIGHT_ON,     
01675          int FOG_ON,           
01676          int SHADOWS_ON,       
01677          int AO_ON,            
01678          int OUTLINE_ON,       
01679          int REFLECTION_ON,    
01680          int TRANSMISSION_ON>  
01681 static __device__ void shader_template(float3 p_obj_color, float3 N, 
01682                                        float p_Ka, float p_Kd, float p_Ks,
01683                                        float p_phong_exp, float p_reflectivity,
01684                                        float p_opacity,
01685                                        float p_outline, float p_outlinewidth,
01686                                        int p_transmode) {
01687   float3 hit_point = ray.origin + t_hit * ray.direction;
01688   float3 result = make_float3(0.0f);
01689   float3 phongcol = make_float3(0.0f);
01690 
01691   // add depth cueing / fog if enabled
01692   // use fog coordinate to modulate importance for AO rays, etc.
01693   float fogmod = 1.0f;
01694   if (FOG_ON && fog_mode != 0) {
01695     fogmod = fog_coord(hit_point);
01696   }
01697 
01698 #if 1
01699   // XXX we really shouldn't have to do this, but it fixes shading 
01700   //     on really bad geometry that can arise from marching cubes 
01701   //     extraction on very noisy cryo-EM and cryo-ET maps
01702   float Ntest = N.x + N.y + N.z;
01703   if (isnan(Ntest) || isinf(Ntest)) {
01704     // add depth cueing / fog if enabled
01705     if (FOG_ON && fogmod < 1.0f) {
01706       result = fog_color(fogmod, result);
01707     }
01708     return;
01709   }
01710 #endif
01711 
01712 #if 1
01713   // don't render transparent surfaces if we've reached the max count
01714   // this implements the same logic as the -trans_max_surfaces argument
01715   // in the CPU version of Tachyon.
01716   if ((p_opacity < 1.0f) && (prd.transcnt < 1)) {
01717     // shoot a transmission ray
01718     PerRayData_radiance new_prd;
01719     new_prd.importance = prd.importance * (1.0f - p_opacity);
01720     new_prd.alpha = 1.0f;
01721     new_prd.result = scene_bg_color;
01722 
01723     // For correct operation with the RTX runtime strategy and its
01724     // associated stack management scheme, we MUST increment the 
01725     // ray recursion depth counter when performing transparent surface
01726     // peeling, otherwise we could go beyond the max recursion depth
01727     // that we previously requested from OptiX.  This will work less well 
01728     // than the former approach in terms of visual outcomes, but we presently
01729     // have no alternative and must avoid issues with runtime stack overruns.
01730     new_prd.depth = prd.depth + 1;
01731     // new_prd.transcnt = max(1, prd.transcnt) - 1; // prevent uint wraparound
01732 
01733     new_prd.transcnt = 0; // don't decrement further since unsigned int type
01734     if (new_prd.importance >= 0.001f && new_prd.depth < max_depth) {
01735 #if defined(ORT_USE_RAY_STEP)
01736 #if defined(ORT_TRANS_USE_INCIDENT)
01737       // step the ray in the incident ray direction
01738       Ray trans_ray = make_Ray(hit_point + ORT_RAY_STEP2, ray.direction, radiance_ray_type, scene_epsilon, RT_DEFAULT_MAX);
01739 #else
01740       // step the ray in the direction opposite the surface normal (going in)
01741       // rather than out, for transmission rays...
01742       Ray trans_ray = make_Ray(hit_point - ORT_RAY_STEP, ray.direction, radiance_ray_type, scene_epsilon, RT_DEFAULT_MAX);
01743 #endif
01744 #else
01745       Ray trans_ray = make_Ray(hit_point, ray.direction, radiance_ray_type, scene_epsilon, RT_DEFAULT_MAX);
01746 #endif
01747       rtTrace(root_object, trans_ray, new_prd);
01748 #if defined(ORT_RAYSTATS)
01749       raystats2_buffer[launch_index].x++; // increment trans ray counter
01750 #endif
01751     }
01752     prd.result = new_prd.result;
01753     return; // early-exit
01754   }
01755 #endif
01756 
01757   // execute the object's texture function
01758   float3 col = p_obj_color; // XXX no texturing implemented yet
01759 
01760   // compute lighting from directional lights
01761 #if defined(VMDOPTIX_LIGHTUSEROBJS)
01762   unsigned int num_dir_lights = dir_light_list.num_lights;
01763   for (int i = 0; i < num_dir_lights; ++i) {
01764     float3 L = dir_light_list.dirs[i];
01765 #else
01766   unsigned int num_dir_lights = dir_lights.size();
01767   for (int i = 0; i < num_dir_lights; ++i) {
01768     DirectionalLight light = dir_lights[i];
01769     float3 L = light.dir;
01770 #endif
01771     shade_light<SHADOWS_ON>(result, hit_point, N, L, p_Kd, p_Ks, p_phong_exp, 
01772                             col, phongcol, RT_DEFAULT_MAX);
01773   }
01774 
01775 #if 0
01776   // compute lighting from positional lights
01777 #if defined(VMDOPTIX_LIGHTUSEROBJS)
01778   unsigned int num_pos_lights = pos_light_list.num_lights;
01779   for (int i = 0; i < num_pos_lights; ++i) {
01780     float3 L = pos_light_list.posns[i];
01781 #else
01782   unsigned int num_pos_lights = pos_lights.size();
01783   for (int i = 0; i < num_pos_lights; ++i) {
01784     PositionalLight light = pos_lights[i];
01785     float3 L = light.pos - hit_point;
01786 #endif
01787     float shadow_tmax = length(L); // compute positional light shadow tmax
01788     L = normalize(L);
01789     shade_light<SHADOWS_ON>(result, hit_point, N, L, p_Kd, p_Ks, p_phong_exp, 
01790                             col, phongcol, shadow_tmax);
01791   }
01792 #endif
01793 
01794   // add point light for camera headlight need for Oculus Rift HMDs,
01795   // equirectangular panorama images, and planetarium dome master images
01796   if (HEADLIGHT_ON && (headlight_mode != 0)) {
01797     float3 L = cam_pos - hit_point;
01798     float shadow_tmax = length(L); // compute positional light shadow tmax
01799     L = normalize(L);
01800     shade_light<SHADOWS_ON>(result, hit_point, N, L, p_Kd, p_Ks, p_phong_exp, 
01801                             col, phongcol, shadow_tmax);
01802   }
01803 
01804   // add ambient occlusion diffuse lighting, if enabled
01805   if (AO_ON && ao_samples > 0) {
01806     result *= ao_direct;
01807     result += ao_ambient * col * p_Kd * shade_ambient_occlusion(hit_point, N, fogmod * p_opacity);
01808   }
01809 
01810   // add edge shading if applicable
01811   if (OUTLINE_ON && p_outline > 0.0f) {
01812     float edgefactor = dot(N, ray.direction);
01813     edgefactor *= edgefactor;
01814     edgefactor = 1.0f - edgefactor;
01815     edgefactor = 1.0f - powf(edgefactor, (1.0f - p_outlinewidth) * 32.0f);
01816     float outlinefactor = __saturatef((1.0f - p_outline) + (edgefactor * p_outline));
01817     result *= outlinefactor;
01818   }
01819 
01820   result += make_float3(p_Ka); // white ambient contribution
01821   result += phongcol;          // add phong highlights
01822 
01823   // spawn reflection rays if necessary
01824   if (REFLECTION_ON && p_reflectivity > 0.0f) {
01825     // ray tree attenuation
01826     PerRayData_radiance new_prd;
01827     new_prd.importance = prd.importance * p_reflectivity;
01828     new_prd.depth = prd.depth + 1;
01829     new_prd.transcnt = prd.transcnt;
01830 
01831     // shoot a reflection ray
01832     if (new_prd.importance >= 0.001f && new_prd.depth < max_depth) {
01833       float3 refl_dir = reflect(ray.direction, N);
01834 #if defined(ORT_USE_RAY_STEP)
01835       Ray refl_ray = make_Ray(hit_point + ORT_RAY_STEP, refl_dir, radiance_ray_type, scene_epsilon, RT_DEFAULT_MAX);
01836 #else
01837       Ray refl_ray = make_Ray(hit_point, refl_dir, radiance_ray_type, scene_epsilon, RT_DEFAULT_MAX);
01838 #endif
01839       rtTrace(root_object, refl_ray, new_prd);
01840 #if defined(ORT_RAYSTATS)
01841       raystats2_buffer[launch_index].w++; // increment refl ray counter
01842 #endif
01843       result += p_reflectivity * new_prd.result;
01844     }
01845   }
01846 
01847   // spawn transmission rays if necessary
01848   float alpha = p_opacity;
01849 
01850 #if 1
01851   if (CLIP_VIEW_ON && (clipview_mode == 2))
01852     sphere_fade_and_clip(hit_point, cam_pos, 
01853                          clipview_start, clipview_end, alpha);
01854 #else
01855   if (CLIP_VIEW_ON && (clipview_mode == 2)) {
01856     // draft implementation of a smooth "fade-out-and-clip sphere"  
01857     float fade_start = 1.00f; // onset of fading 
01858     float fade_end   = 0.20f; // fully transparent
01859     float camdist = length(hit_point - cam_pos);
01860 
01861     // XXX we can omit the distance test since alpha modulation value is clamped
01862     // if (1 || camdist < fade_start) {
01863       float fade_len = fade_start - fade_end;
01864       alpha *= __saturatef((camdist - fade_start) / fade_len);
01865     // }
01866   }
01867 #endif
01868 
01869   // TRANSMISSION_ON: handles transparent surface shading, test is only
01870   // performed when the VMD geometry has a known-transparent material
01871   // CLIP_VIEW_ON: forces check of alpha value for all geom as per transparent 
01872   // material, since all geometry may become tranparent with the 
01873   // fade+clip sphere active
01874 #if defined(CADENSMOVIE)
01875   // XXX specially modified shading of transparent surfaces used to 
01876   //     give "cloud like" appearances to molecular surfaces 
01877   if ((p_transmode>1 || (TRANSMISSION_ON || CLIP_VIEW_ON)) && alpha < 0.999f ) {
01878     // Emulate Tachyon/Raster3D's angle-dependent surface opacity if enabled
01879     switch(int(p_transmode)) {
01880 
01881     case 3: // transmode=3 ==> alternate angle-modulated-transparency - try to not let it be fully opaque at edges
01882         alpha = 1.0f + cosf(3.1415926f * (1.0f-alpha) * dot(N, ray.direction));
01883         alpha = alpha*alpha * 0.165;  // something less than 0.25
01884         break;
01885 
01886     case 2: {
01887         // Glow case.  Why can't we make this be case 2, and set our special material's transmode to 2??
01888         float d = fminf( dot(N, ray.direction), 0.5f );  // flat-topped glow
01889         float v = alpha * fabsf(d*d); // - 0.5f;
01890         result = 3.5f * v * result; // * make_float3( v,v,v );
01891         alpha = v*0.15;
01892         // alpha = v*v*0.25f;  
01893         // result = make_float3( alpha, p_transmode*(0.5f - fabsf(d*d)), 0 );
01894         break;
01895       }
01896 
01897     case 1:
01898         // normal "angle-modulated transparency" case
01899         alpha = 1.0f + cosf(3.1415926f * (1.0f-alpha) * dot(N, ray.direction));
01900         alpha = alpha*alpha * 0.25f;
01901         break;
01902     }
01903 #else
01904   if ((TRANSMISSION_ON || CLIP_VIEW_ON) && alpha < 0.999f ) {
01905     // Emulate Tachyon/Raster3D's angle-dependent surface opacity if enabled
01906     if (p_transmode) {
01907       alpha = 1.0f + cosf(3.1415926f * (1.0f-alpha) * dot(N, ray.direction));
01908       alpha = alpha*alpha * 0.25f;
01909     }
01910 #endif
01911 
01912     result *= alpha; // scale down current lighting by opacity
01913 
01914     // shoot a transmission ray
01915     PerRayData_radiance new_prd;
01916     new_prd.importance = prd.importance * (1.0f - alpha);
01917     new_prd.alpha = 1.0f;
01918     new_prd.result = scene_bg_color;
01919     new_prd.depth = prd.depth + 1;
01920     new_prd.transcnt = max(1, prd.transcnt) - 1; // prevent uint wraparound
01921     if (new_prd.importance >= 0.001f && new_prd.depth < max_depth) {
01922 #if defined(ORT_USE_RAY_STEP)
01923 #if defined(ORT_TRANS_USE_INCIDENT)
01924       // step the ray in the incident ray direction
01925       Ray trans_ray = make_Ray(hit_point + ORT_RAY_STEP2, ray.direction, radiance_ray_type, scene_epsilon, RT_DEFAULT_MAX);
01926 #else
01927       // step the ray in the direction opposite the surface normal (going in)
01928       // rather than out, for transmission rays...
01929       Ray trans_ray = make_Ray(hit_point - ORT_RAY_STEP, ray.direction, radiance_ray_type, scene_epsilon, RT_DEFAULT_MAX);
01930 #endif
01931 #else
01932       Ray trans_ray = make_Ray(hit_point, ray.direction, radiance_ray_type, scene_epsilon, RT_DEFAULT_MAX);
01933 #endif
01934       rtTrace(root_object, trans_ray, new_prd);
01935 #if defined(ORT_RAYSTATS)
01936       raystats2_buffer[launch_index].x++; // increment trans ray counter
01937 #endif
01938     }
01939     result += (1.0f - alpha) * new_prd.result;
01940     prd.alpha = alpha + (1.0f - alpha) * new_prd.alpha; 
01941   }
01942 
01943   // add depth cueing / fog if enabled
01944   if (FOG_ON && fogmod < 1.0f) {
01945     result = fog_color(fogmod, result);
01946   }
01947 
01948   prd.result = result; // pass the color back up the tree
01949 }
01950 
01951 
01952 // color state associated with meshes or other primitives that 
01953 // don't provide per-vertex, or per-facet color data
01954 rtDeclareVariable(float3, uniform_color, , );
01955 
01956 // VMD material shading coefficients
01957 rtDeclareVariable(float, Ka, , );
01958 rtDeclareVariable(float, Kd, , );
01959 rtDeclareVariable(float, Ks, , );
01960 rtDeclareVariable(float, phong_exp, , );
01961 rtDeclareVariable(float, Krefl, , );
01962 rtDeclareVariable(float, opacity, , );
01963 rtDeclareVariable(float, outline, , );
01964 rtDeclareVariable(float, outlinewidth, , );
01965 rtDeclareVariable(int, transmode, , );
01966 
01967 
01968 // Any hit program for opaque objects
01969 RT_PROGRAM void any_hit_opaque() {
01970   // this material is opaque, so it fully attenuates all shadow rays
01971   prd_shadow.attenuation = optix::make_float3(0.0f);
01972 
01973   // upon hitting _any_ opaque object during an anyhit traversal, 
01974   // we should immediately terminate traversal
01975   rtTerminateRay();
01976 }
01977 
01978 
01979 // Any hit program required for shadow filtering through transparent materials
01980 RT_PROGRAM void any_hit_transmission() {
01981   // use a VERY simple shadow filtering scheme based on opacity
01982   prd_shadow.attenuation *= make_float3(1.0f - opacity);
01983 
01984   // check to see if we've hit 100% shadow or not
01985   if (fmaxf(prd_shadow.attenuation) < 0.001f) {
01986     rtTerminateRay();
01987   } else {
01988 #if defined(ORT_RAYSTATS)
01989     raystats2_buffer[launch_index].y++; // increment trans ray skip counter
01990 #endif
01991     rtIgnoreIntersection();
01992   }
01993 }
01994 
01995 
01996 // Any hit program required for shadow filtering when an 
01997 // HMD/camera fade-and-clip is active, through both 
01998 // solid and transparent materials
01999 RT_PROGRAM void any_hit_clip_sphere() {
02000 
02001   // compute hit point for use in evaluating fade/clip effect
02002   float3 hit_point = ray.origin + t_hit * ray.direction;
02003 
02004   // compute additional attenuation from clipping sphere if enabled
02005   float clipalpha = 1.0f;
02006   if (clipview_mode == 2) {
02007     sphere_fade_and_clip(hit_point, cam_pos, clipview_start, clipview_end, 
02008                          clipalpha);
02009   }
02010 
02011   // use a VERY simple shadow filtering scheme based on opacity
02012   prd_shadow.attenuation = make_float3(1.0f - (clipalpha * opacity));
02013 
02014   // check to see if we've hit 100% shadow or not
02015   if (fmaxf(prd_shadow.attenuation) < 0.001f) {
02016     rtTerminateRay();
02017   } else {
02018 #if defined(ORT_RAYSTATS)
02019     raystats2_buffer[launch_index].y++; // increment trans ray skip counter
02020 #endif
02021     rtIgnoreIntersection();
02022   }
02023 }
02024 
02025 
02026 
02027 // normal calc routine needed only to simplify the macro to produce the
02028 // complete combinatorial expansion of template-specialized 
02029 // closest hit radiance functions 
02030 static __inline__ __device__ float3 calc_ffworld_normal(const float3 &Nshading, const float3 &Ngeometric) {
02031   float3 world_shading_normal = normalize(rtTransformNormal(RT_OBJECT_TO_WORLD, Nshading));
02032   float3 world_geometric_normal = normalize(rtTransformNormal(RT_OBJECT_TO_WORLD, Ngeometric));
02033   return faceforward(world_shading_normal, -ray.direction, world_geometric_normal);
02034 }
02035 
02036 
02037 template<int HWTRIANGLES> static __device__ __inline__ 
02038 float3 get_intersection_normal() {
02039 #if defined(ORT_USERTXAPIS)
02040   // OptiX RTX hardware triangle API
02041   if (HWTRIANGLES) {
02042     const unsigned int primIdx = rtGetPrimitiveIndex();
02043     const float3 Ng = unpackNormal(normalBuffer[primIdx].x);
02044     float3 Ns;
02045     if (has_vertex_normals) {
02046       const float3& n0 = unpackNormal(normalBuffer[primIdx].y);
02047       const float3& n1 = unpackNormal(normalBuffer[primIdx].z);
02048       const float3& n2 = unpackNormal(normalBuffer[primIdx].w);
02049       Ns = optix::normalize(n0 * (1.0f - barycentrics.x - barycentrics.y) +
02050                             n1 * barycentrics.x + 
02051                             n2 * barycentrics.y);
02052     } else {
02053       Ns = Ng;
02054     }
02055     return calc_ffworld_normal(Ns, Ng);
02056   } else {
02057     // classic OptiX APIs and non-triangle geometry
02058     return calc_ffworld_normal(shading_normal, geometric_normal);
02059   } 
02060 #else
02061   // classic OptiX APIs and non-triangle geometry
02062   return calc_ffworld_normal(shading_normal, geometric_normal);
02063 #endif
02064 }
02065 
02066 
02067 template<int HWTRIANGLES> static __device__ __inline__ 
02068 float3 get_intersection_color() {
02069 #if defined(ORT_USERTXAPIS)
02070   // OptiX RTX hardware triangle API
02071   if (HWTRIANGLES) {
02072     if (has_vertex_colors) {
02073       const float ci2f = 1.0f / 255.0f;
02074       const unsigned int primIdx = rtGetPrimitiveIndex();
02075       const float3 c0 = colorBuffer[3 * primIdx + 0] * ci2f;
02076       const float3 c1 = colorBuffer[3 * primIdx + 1] * ci2f;
02077       const float3 c2 = colorBuffer[3 * primIdx + 2] * ci2f;
02078 
02079       // interpolate triangle color from barycentrics
02080       return (c0 * (1.0f - barycentrics.x - barycentrics.y) +
02081               c1 * barycentrics.x + 
02082               c2 * barycentrics.y);
02083     } else {
02084       return uniform_color; // return uniform mesh color
02085     }
02086   } else {
02087     // classic OptiX APIs and non-triangle geometry
02088     return obj_color; // return object color determined during intersection
02089   }
02090 #else
02091   // classic OptiX APIs and non-triangle geometry
02092   return obj_color; // return object color determined during intersection
02093 #endif
02094 }
02095 
02096 
02097 // general-purpose any-hit program, with all template options enabled, 
02098 // intended for shader debugging and comparison with the original
02099 // Tachyon full_shade() code.
02100 RT_PROGRAM void closest_hit_radiance_general() {
02101   shader_template<1, 1, 1, 1, 1, 1, 1, 1>(get_intersection_color<0>(),
02102                                           get_intersection_normal<0>(), 
02103                                           Ka, Kd, Ks, phong_exp, Krefl, opacity,
02104                                           outline, outlinewidth, transmode);
02105 }
02106 
02107 #if defined(ORT_USERTXAPIS)
02108 // OptiX RTX hardware triangle API
02109 RT_PROGRAM void closest_hit_radiance_general_hwtri() {
02110   shader_template<1, 1, 1, 1, 1, 1, 1, 1>(get_intersection_color<1>(),
02111                                           get_intersection_normal<1>(), 
02112                                           Ka, Kd, Ks, phong_exp, Krefl, opacity,
02113                                           outline, outlinewidth, transmode);
02114 }
02115 #endif
02116 
02117 
02118 //
02119 // Object and/or vertex/color/normal buffers...
02120 //
02121 
02122 // cylinder array buffers
02123 rtBuffer<vmd_cylinder> cylinder_buffer;
02124 rtBuffer<vmd_cylinder_color> cylinder_color_buffer;
02125 
02126 // ring array buffer
02127 rtBuffer<vmd_ring_color> ring_color_buffer;
02128 
02129 // sphere array buffers
02130 rtBuffer<vmd_sphere> sphere_buffer;
02131 rtBuffer<vmd_sphere_color> sphere_color_buffer;
02132 
02133 // triangle mesh buffers
02134 rtBuffer<vmd_tricolor> tricolor_buffer;
02135 rtBuffer<vmd_trimesh_c4u_n3b_v3f> trimesh_c4u_n3b_v3f_buffer;
02136 rtBuffer<vmd_trimesh_n3f_v3f> trimesh_n3f_v3f_buffer;
02137 rtBuffer<vmd_trimesh_n3b_v3f> trimesh_n3b_v3f_buffer;
02138 rtBuffer<vmd_trimesh_v3f> trimesh_v3f_buffer;
02139 
02140 
02141 //
02142 // Cylinder array primitive
02143 //
02144 RT_PROGRAM void cylinder_array_intersect(int primIdx) {
02145   float3 start = cylinder_buffer[primIdx].start;
02146   float radius = cylinder_buffer[primIdx].radius;
02147   float3 axis = cylinder_buffer[primIdx].axis;
02148 
02149   float3 rc = ray.origin - start;
02150   float3 n = cross(ray.direction, axis);
02151   float ln = length(n);
02152 
02153   // check if ray is parallel to cylinder
02154   if (ln == 0.0f) {
02155     return; // ray is parallel, we missed or went through the "hole"
02156   } 
02157   n /= ln;
02158   float d = fabsf(dot(rc, n));
02159 
02160   // check for cylinder intersection
02161   if (d <= radius) {
02162     float3 O = cross(rc, axis);
02163     float t = -dot(O, n) / ln;
02164     O = cross(n, axis);
02165     O = normalize(O);
02166     float s = dot(ray.direction, O); 
02167     s = fabs(sqrtf(radius*radius - d*d) / s);
02168     float axlen = length(axis);
02169     float3 axis_u = normalize(axis);
02170 
02171     // test hit point against cylinder ends
02172     float tin = t - s;
02173     float3 hit = ray.origin + ray.direction * tin;
02174     float3 tmp2 = hit - start;
02175     float tmp = dot(tmp2, axis_u);
02176     if ((tmp > 0.0f) && (tmp < axlen)) {
02177       if (rtPotentialIntersection(tin)) {
02178         shading_normal = geometric_normal = normalize(hit - (tmp * axis_u + start));
02179 
02180         // uniform color for the entire object
02181         obj_color = uniform_color;
02182         rtReportIntersection(0);
02183       }
02184     }
02185     
02186     // continue with second test...
02187     float tout = t + s;
02188     hit = ray.origin + ray.direction * tout;
02189     tmp2 = hit - start;
02190     tmp = dot(tmp2, axis_u);
02191     if ((tmp > 0.0f) && (tmp < axlen)) {
02192       if (rtPotentialIntersection(tout)) {
02193         shading_normal = geometric_normal = normalize(hit - (tmp * axis_u + start));
02194 
02195         // uniform color for the entire object
02196         obj_color = uniform_color;
02197         rtReportIntersection(0);
02198       }
02199     }
02200   }
02201 }
02202 
02203 
02204 RT_PROGRAM void cylinder_array_bounds(int primIdx, float result[6]) {
02205   const float3 start = cylinder_buffer[primIdx].start;
02206   const float3 end = start + cylinder_buffer[primIdx].axis;
02207   const float3 rad = make_float3(cylinder_buffer[primIdx].radius);
02208   optix::Aabb* aabb = (optix::Aabb*)result;
02209 
02210   if (rad.x > 0.0f && !isinf(rad.x)) {
02211     aabb->m_min = fminf(start - rad, end - rad);
02212     aabb->m_max = fmaxf(start + rad, end + rad);
02213   } else {
02214     aabb->invalidate();
02215   }
02216 }
02217 
02218 
02219 //
02220 // Color-per-cylinder array primitive
02221 //
02222 RT_PROGRAM void cylinder_array_color_intersect(int primIdx) {
02223   float3 start = cylinder_color_buffer[primIdx].start;
02224   float radius = cylinder_color_buffer[primIdx].radius;
02225   float3 axis = cylinder_color_buffer[primIdx].axis;
02226 
02227   float3 rc = ray.origin - start;
02228   float3 n = cross(ray.direction, axis);
02229   float ln = length(n);
02230 
02231   // check if ray is parallel to cylinder
02232   if (ln == 0.0f) {
02233     return; // ray is parallel, we missed or went through the "hole"
02234   } 
02235   n /= ln;
02236   float d = fabsf(dot(rc, n));
02237 
02238   // check for cylinder intersection
02239   if (d <= radius) {
02240     float3 O = cross(rc, axis);
02241     float t = -dot(O, n) / ln;
02242     O = cross(n, axis);
02243     O = normalize(O);
02244     float s = dot(ray.direction, O); 
02245     s = fabs(sqrtf(radius*radius - d*d) / s);
02246     float axlen = length(axis);
02247     float3 axis_u = normalize(axis);
02248 
02249     // test hit point against cylinder ends
02250     float tin = t - s;
02251     float3 hit = ray.origin + ray.direction * tin;
02252     float3 tmp2 = hit - start;
02253     float tmp = dot(tmp2, axis_u);
02254     if ((tmp > 0.0f) && (tmp < axlen)) {
02255       if (rtPotentialIntersection(tin)) {
02256         shading_normal = geometric_normal = normalize(hit - (tmp * axis_u + start));
02257         obj_color = cylinder_color_buffer[primIdx].color;
02258         rtReportIntersection(0);
02259       }
02260     }
02261     
02262     // continue with second test...
02263     float tout = t + s;
02264     hit = ray.origin + ray.direction * tout;
02265     tmp2 = hit - start;
02266     tmp = dot(tmp2, axis_u);
02267     if ((tmp > 0.0f) && (tmp < axlen)) {
02268       if (rtPotentialIntersection(tout)) {
02269         shading_normal = geometric_normal = normalize(hit - (tmp * axis_u + start));
02270         obj_color = cylinder_color_buffer[primIdx].color;
02271         rtReportIntersection(0);
02272       }
02273     }
02274   }
02275 }
02276 
02277 
02278 RT_PROGRAM void cylinder_array_color_bounds(int primIdx, float result[6]) {
02279   const float3 start = cylinder_color_buffer[primIdx].start;
02280   const float3 end = start + cylinder_color_buffer[primIdx].axis;
02281   const float3 rad = make_float3(cylinder_color_buffer[primIdx].radius);
02282   optix::Aabb* aabb = (optix::Aabb*)result;
02283 
02284   if (rad.x > 0.0f && !isinf(rad.x)) {
02285     aabb->m_min = fminf(start - rad, end - rad);
02286     aabb->m_max = fmaxf(start + rad, end + rad);
02287   } else {
02288     aabb->invalidate();
02289   }
02290 }
02291 
02292 
02293 //
02294 // Ring array primitive
02295 //
02296 RT_PROGRAM void ring_array_color_intersect(int primIdx) {
02297   const float3 center = ring_color_buffer[primIdx].center;
02298   const float3 norm = ring_color_buffer[primIdx].norm;
02299   const float inrad = ring_color_buffer[primIdx].inrad;
02300   const float outrad = ring_color_buffer[primIdx].outrad;
02301   const float3 color = ring_color_buffer[primIdx].color;
02302 
02303   float d = -dot(center, norm); 
02304   float t = -(d + dot(norm, ray.origin));
02305   float td = dot(norm, ray.direction);
02306   if (td != 0.0f) {
02307     t /= td;
02308     if (t >= 0.0f) {
02309       float3 hit = ray.origin + t * ray.direction;
02310       float rd = length(hit - center);
02311       if ((rd > inrad) && (rd < outrad)) {
02312         if (rtPotentialIntersection(t)) {
02313           shading_normal = geometric_normal = norm;
02314           obj_color = color;
02315           rtReportIntersection(0);
02316         }
02317       }
02318     }
02319   }
02320 }
02321 
02322 
02323 RT_PROGRAM void ring_array_color_bounds(int primIdx, float result[6]) {
02324   const float3 center = ring_color_buffer[primIdx].center;
02325   const float3 rad = make_float3(ring_color_buffer[primIdx].outrad);
02326   optix::Aabb* aabb = (optix::Aabb*)result;
02327 
02328   if (rad.x > 0.0f && !isinf(rad.x)) {
02329     aabb->m_min = center - rad;
02330     aabb->m_max = center + rad;
02331   } else {
02332     aabb->invalidate();
02333   }
02334 }
02335 
02336 
02337 
02338 #if defined(ORT_USE_SPHERES_HEARNBAKER)
02339 
02340 // Ray-sphere intersection method with improved floating point precision
02341 // for cases where the sphere size is small relative to the distance
02342 // from the camera to the sphere.  This implementation is based on
02343 // Eq. 10-72, p.603 of "Computer Graphics with OpenGL", 3rd Ed.,
02344 // by Donald Hearn and Pauline Baker, 2004.  Shown in Eq. 10, p.639
02345 // in the 4th edition of the book (Hearn, Baker, Carithers).
02346 static __host__ __device__ __inline__
02347 void sphere_intersect_hearn_baker2(float3 center, float radius,
02348                                    const float3 &spcolor) {
02349   float3 deltap = center - ray.origin;
02350   float ddp = dot(ray.direction, deltap);
02351   float3 remedyTerm = deltap - ddp * ray.direction;
02352   float disc = radius*radius - dot(remedyTerm, remedyTerm);
02353   if (disc >= 0.0f) {
02354     float disc_root = sqrtf(disc);
02355     float t1 = ddp - disc_root;
02356     float t2 = ddp + disc_root;
02357 
02358     if (rtPotentialIntersection(t1)) {
02359       shading_normal = geometric_normal = (t1*ray.direction - deltap) / radius;
02360       obj_color = spcolor;
02361       rtReportIntersection(0);
02362     }
02363 
02364     if (rtPotentialIntersection(t2)) {
02365       shading_normal = geometric_normal = (t2*ray.direction - deltap) / radius;
02366       obj_color = spcolor;
02367       rtReportIntersection(0);
02368     }
02369   }
02370 }
02371 
02372 #else
02373 
02374 //
02375 // Ray-sphere intersection using standard geometric solution approach
02376 //
02377 static __host__ __device__ __inline__
02378 void sphere_intersect_classic(float3 center, float radius,
02379                               const float3 &spcolor) {
02380   float3 V = center - ray.origin;
02381   float b = dot(V, ray.direction);
02382   float disc = b*b + radius*radius - dot(V, V);
02383   if (disc > 0.0f) {
02384     disc = sqrtf(disc);
02385 
02386 //#define FASTONESIDEDSPHERES 1
02387 #if defined(FASTONESIDEDSPHERES)
02388     // only calculate the nearest intersection, for speed
02389     float t1 = b - disc;
02390     if (rtPotentialIntersection(t1)) {
02391       shading_normal = geometric_normal = (t1*ray.direction - V) / radius;
02392       obj_color = spcolor;
02393       rtReportIntersection(0);
02394     }
02395 #else
02396     float t2 = b + disc;
02397     if (rtPotentialIntersection(t2)) {
02398       shading_normal = geometric_normal = (t2*ray.direction - V) / radius;
02399       obj_color = spcolor;
02400       rtReportIntersection(0);
02401     }
02402 
02403     float t1 = b - disc;
02404     if (rtPotentialIntersection(t1)) {
02405       shading_normal = geometric_normal = (t1*ray.direction - V) / radius;
02406       obj_color = spcolor;
02407       rtReportIntersection(0);
02408     }
02409 #endif
02410   }
02411 }
02412 
02413 #endif
02414 
02415 
02416 //
02417 // Sphere array primitive
02418 //
02419 RT_PROGRAM void sphere_array_intersect(int primIdx) {
02420   float3 center = sphere_buffer[primIdx].center;
02421   float radius = sphere_buffer[primIdx].radius;
02422 
02423   // uniform color for the entire object
02424 #if defined(ORT_USE_SPHERES_HEARNBAKER)
02425   sphere_intersect_hearn_baker2(center, radius, uniform_color);
02426 #else
02427   sphere_intersect_classic(center, radius, uniform_color);
02428 #endif
02429 }
02430 
02431 
02432 RT_PROGRAM void sphere_array_bounds(int primIdx, float result[6]) {
02433   const float3 cen = sphere_buffer[primIdx].center;
02434   const float3 rad = make_float3(sphere_buffer[primIdx].radius);
02435   optix::Aabb* aabb = (optix::Aabb*)result;
02436 
02437   if (rad.x > 0.0f && !isinf(rad.x)) {
02438     aabb->m_min = cen - rad;
02439     aabb->m_max = cen + rad;
02440   } else {
02441     aabb->invalidate();
02442   }
02443 }
02444 
02445 
02446 //
02447 // Color-per-sphere sphere array 
02448 //
02449 RT_PROGRAM void sphere_array_color_intersect(int primIdx) {
02450   float3 center = sphere_color_buffer[primIdx].center;
02451   float radius = sphere_color_buffer[primIdx].radius;
02452 
02453 #if defined(ORT_USE_SPHERES_HEARNBAKER)
02454   sphere_intersect_hearn_baker2(center, radius, sphere_color_buffer[primIdx].color);
02455 #else
02456   sphere_intersect_classic(center, radius, sphere_color_buffer[primIdx].color);
02457 #endif
02458 }
02459 
02460 
02461 RT_PROGRAM void sphere_array_color_bounds(int primIdx, float result[6]) {
02462   const float3 cen = sphere_color_buffer[primIdx].center;
02463   const float3 rad = make_float3(sphere_color_buffer[primIdx].radius);
02464   optix::Aabb* aabb = (optix::Aabb*)result;
02465 
02466   if (rad.x > 0.0f && !isinf(rad.x)) {
02467     aabb->m_min = cen - rad;
02468     aabb->m_max = cen + rad;
02469   } else {
02470     aabb->invalidate();
02471   }
02472 }
02473 
02474 
02475 
02476 //
02477 // Triangle list primitive - unstructured triangle soup
02478 //
02479 
02480 
02481 // inline device function for computing triangle bounding boxes
02482 __device__ __inline__ void generic_tri_bounds(optix::Aabb *aabb,
02483                                               float3 v0, float3 v1, float3 v2) {
02484 #if 1
02485   // conventional paranoid implementation that culls degenerate triangles
02486   float area = length(cross(v1-v0, v2-v0));
02487   if (area > 0.0f && !isinf(area)) {
02488     aabb->m_min = fminf(fminf(v0, v1), v2);
02489     aabb->m_max = fmaxf(fmaxf(v0, v1), v2);
02490   } else {
02491     aabb->invalidate();
02492   }
02493 #else
02494   // don't cull any triangles, even if they might be degenerate
02495   aabb->m_min = fminf(fminf(v0, v1), v2);
02496   aabb->m_max = fmaxf(fmaxf(v0, v1), v2);
02497 #endif
02498 }
02499 
02500 
02501 
02502 RT_PROGRAM void tricolor_intersect(int primIdx) {
02503   float3 v0 = tricolor_buffer[primIdx].v0;
02504   float3 v1 = tricolor_buffer[primIdx].v1;
02505   float3 v2 = tricolor_buffer[primIdx].v2;
02506 
02507   // Intersect ray with triangle
02508   float3 n;
02509   float t, beta, gamma;
02510   if (intersect_triangle(ray, v0, v1, v2, n, t, beta, gamma)) {
02511     if (rtPotentialIntersection(t)) {
02512       float3 n0 = tricolor_buffer[primIdx].n0;
02513       float3 n1 = tricolor_buffer[primIdx].n1;
02514       float3 n2 = tricolor_buffer[primIdx].n2;
02515       shading_normal = normalize(n1*beta + n2*gamma + n0*(1.0f-beta-gamma));
02516       geometric_normal = normalize(n);
02517 
02518       float3 c0 = tricolor_buffer[primIdx].c0;
02519       float3 c1 = tricolor_buffer[primIdx].c1;
02520       float3 c2 = tricolor_buffer[primIdx].c2;
02521       obj_color = c1*beta + c2*gamma + c0*(1.0f-beta-gamma);
02522       rtReportIntersection(0);
02523     }
02524   }
02525 }
02526 
02527 RT_PROGRAM void tricolor_bounds(int primIdx, float result[6]) {
02528   float3 v0 = tricolor_buffer[primIdx].v0;
02529   float3 v1 = tricolor_buffer[primIdx].v1;
02530   float3 v2 = tricolor_buffer[primIdx].v2;
02531 
02532   optix::Aabb *aabb = (optix::Aabb*)result;
02533   generic_tri_bounds(aabb, v0, v1, v2);
02534 }
02535 
02536 
02537 
02538 RT_PROGRAM void trimesh_c4u_n3b_v3f_intersect(int primIdx) {
02539   float3 v0 = trimesh_c4u_n3b_v3f_buffer[primIdx].v0;
02540   float3 v1 = trimesh_c4u_n3b_v3f_buffer[primIdx].v1;
02541   float3 v2 = trimesh_c4u_n3b_v3f_buffer[primIdx].v2;
02542 
02543   // Intersect ray with triangle
02544   float3 n;
02545   float t, beta, gamma;
02546   if (intersect_triangle(ray, v0, v1, v2, n, t, beta, gamma)) {
02547     if (rtPotentialIntersection(t)) {
02548       // conversion from GLbyte format, Table 2.6, p. 44 of OpenGL spec 1.2.1
02549       // float = (2c+1)/(2^8-1)
02550       const float ci2f = 1.0f / 255.0f;
02551       const float cn2f = 1.0f / 127.5f;
02552 
02553       float3 n0 = trimesh_c4u_n3b_v3f_buffer[primIdx].n0 * cn2f + ci2f;
02554       float3 n1 = trimesh_c4u_n3b_v3f_buffer[primIdx].n1 * cn2f + ci2f;
02555       float3 n2 = trimesh_c4u_n3b_v3f_buffer[primIdx].n2 * cn2f + ci2f;
02556       shading_normal = normalize(n1*beta + n2*gamma + n0*(1.0f-beta-gamma));
02557       geometric_normal = normalize(n);
02558 
02559       // conversion from GLubyte format, Table 2.6, p. 44 of OpenGL spec 1.2.1
02560       // float = c/(2^8-1)
02561       float3 c0 = trimesh_c4u_n3b_v3f_buffer[primIdx].c0 * ci2f;
02562       float3 c1 = trimesh_c4u_n3b_v3f_buffer[primIdx].c1 * ci2f;
02563       float3 c2 = trimesh_c4u_n3b_v3f_buffer[primIdx].c2 * ci2f;
02564       obj_color = c1*beta + c2*gamma + c0*(1.0f-beta-gamma);
02565       rtReportIntersection(0);
02566     }
02567   }
02568 }
02569 
02570 RT_PROGRAM void trimesh_c4u_n3b_v3f_bounds(int primIdx, float result[6]) {
02571   float3 v0 = trimesh_c4u_n3b_v3f_buffer[primIdx].v0;
02572   float3 v1 = trimesh_c4u_n3b_v3f_buffer[primIdx].v1;
02573   float3 v2 = trimesh_c4u_n3b_v3f_buffer[primIdx].v2;
02574 
02575   optix::Aabb *aabb = (optix::Aabb*)result;
02576   generic_tri_bounds(aabb, v0, v1, v2);
02577 }
02578 
02579 
02580 
02581 RT_PROGRAM void trimesh_n3f_v3f_intersect(int primIdx) {
02582   float3 v0 = trimesh_n3f_v3f_buffer[primIdx].v0;
02583   float3 v1 = trimesh_n3f_v3f_buffer[primIdx].v1;
02584   float3 v2 = trimesh_n3f_v3f_buffer[primIdx].v2;
02585 
02586   // Intersect ray with triangle
02587   float3 n;
02588   float t, beta, gamma;
02589   if (intersect_triangle(ray, v0, v1, v2, n, t, beta, gamma)) {
02590     if (rtPotentialIntersection(t)) {
02591       float3 n0 = trimesh_n3f_v3f_buffer[primIdx].n0;
02592       float3 n1 = trimesh_n3f_v3f_buffer[primIdx].n1;
02593       float3 n2 = trimesh_n3f_v3f_buffer[primIdx].n2;
02594 
02595       shading_normal = normalize(n1*beta + n2*gamma + n0*(1.0f-beta-gamma) );
02596       geometric_normal = normalize(n);
02597 
02598       // uniform color for the entire object
02599       obj_color = uniform_color;
02600       rtReportIntersection(0);
02601     }
02602   }
02603 }
02604 
02605 RT_PROGRAM void trimesh_n3f_v3f_bounds(int primIdx, float result[6]) {
02606   float3 v0 = trimesh_n3f_v3f_buffer[primIdx].v0;
02607   float3 v1 = trimesh_n3f_v3f_buffer[primIdx].v1;
02608   float3 v2 = trimesh_n3f_v3f_buffer[primIdx].v2;
02609 
02610   optix::Aabb *aabb = (optix::Aabb*)result;
02611   generic_tri_bounds(aabb, v0, v1, v2);
02612 }
02613 
02614 
02615 
02616 RT_PROGRAM void trimesh_n3b_v3f_intersect(int primIdx) {
02617   float3 v0 = trimesh_n3b_v3f_buffer[primIdx].v0;
02618   float3 v1 = trimesh_n3b_v3f_buffer[primIdx].v1;
02619   float3 v2 = trimesh_n3b_v3f_buffer[primIdx].v2;
02620 
02621   // Intersect ray with triangle
02622   float3 n;
02623   float t, beta, gamma;
02624   if (intersect_triangle(ray, v0, v1, v2, n, t, beta, gamma)) {
02625     if (rtPotentialIntersection(t)) {
02626       // conversion from GLbyte format, Table 2.6, p. 44 of OpenGL spec 1.2.1
02627       // float = (2c+1)/(2^8-1)
02628       const float ci2f = 1.0f / 255.0f;
02629       const float cn2f = 1.0f / 127.5f;
02630 
02631       float3 n0 = trimesh_n3b_v3f_buffer[primIdx].n0 * cn2f + ci2f;
02632       float3 n1 = trimesh_n3b_v3f_buffer[primIdx].n1 * cn2f + ci2f;
02633       float3 n2 = trimesh_n3b_v3f_buffer[primIdx].n2 * cn2f + ci2f;
02634       shading_normal = normalize(n1*beta + n2*gamma + n0*(1.0f-beta-gamma) );
02635       geometric_normal = normalize(n);
02636 
02637       // uniform color for the entire object
02638       obj_color = uniform_color;
02639       rtReportIntersection(0);
02640     }
02641   }
02642 }
02643 
02644 RT_PROGRAM void trimesh_n3b_v3f_bounds(int primIdx, float result[6]) {
02645   float3 v0 = trimesh_n3b_v3f_buffer[primIdx].v0;
02646   float3 v1 = trimesh_n3b_v3f_buffer[primIdx].v1;
02647   float3 v2 = trimesh_n3b_v3f_buffer[primIdx].v2;
02648 
02649   optix::Aabb *aabb = (optix::Aabb*)result;
02650   generic_tri_bounds(aabb, v0, v1, v2);
02651 }
02652 
02653 
02654 RT_PROGRAM void trimesh_v3f_intersect(int primIdx) {
02655   float3 v0 = trimesh_v3f_buffer[primIdx].v0;
02656   float3 v1 = trimesh_v3f_buffer[primIdx].v1;
02657   float3 v2 = trimesh_v3f_buffer[primIdx].v2;
02658 
02659   // Intersect ray with triangle
02660   float3 n;
02661   float t, beta, gamma;
02662   if (intersect_triangle(ray, v0, v1, v2, n, t, beta, gamma)) {
02663     if (rtPotentialIntersection(t)) {
02664       shading_normal = geometric_normal = normalize(n);
02665 
02666       // uniform color for the entire object
02667       obj_color = uniform_color;
02668       rtReportIntersection(0);
02669     }
02670   }
02671 }
02672 
02673 RT_PROGRAM void trimesh_v3f_bounds(int primIdx, float result[6]) {
02674   float3 v0 = trimesh_v3f_buffer[primIdx].v0;
02675   float3 v1 = trimesh_v3f_buffer[primIdx].v1;
02676   float3 v2 = trimesh_v3f_buffer[primIdx].v2;
02677 
02678   optix::Aabb *aabb = (optix::Aabb*)result;
02679   generic_tri_bounds(aabb, v0, v1, v2);
02680 }
02681 
02682 
02683 //
02684 // The code below this point generates the combinatorial expansion of
02685 // the series of shading feature on/off flags.  This produces a complete
02686 // set of template-specialized shader variants that handle every 
02687 // performance-beneficial combination of scene-wide and material-specific 
02688 // shading features, which are then linked to the OptiX scene graph 
02689 // material nodes.
02690 //
02691 #if defined(ORT_USE_TEMPLATE_SHADERS)
02692 
02693 #define CLIP_VIEW_on      1
02694 #define CLIP_VIEW_off     0
02695 #define HEADLIGHT_on      1
02696 #define HEADLIGHT_off     0
02697 #define FOG_on            1
02698 #define FOG_off           0
02699 #define SHADOWS_on        1
02700 #define SHADOWS_off       0
02701 #define AO_on             1
02702 #define AO_off            0
02703 #define OUTLINE_on        1
02704 #define OUTLINE_off       0
02705 #define REFL_on           1
02706 #define REFL_off          0
02707 #define TRANS_on          1
02708 #define TRANS_off         0
02709 
02710 #define DEFINE_CLOSEST_HIT( mclipview, mhlight, mfog, mshad, mao, moutl, mrefl, mtrans ) \
02711   RT_PROGRAM void                                                        \
02712   closest_hit_radiance_CLIP_VIEW_##mclipview##_HEADLIGHT_##mhlight##_FOG_##mfog##_SHADOWS_##mshad##_AO_##mao##_OUTLINE_##moutl##_REFL_##mrefl##_TRANS_##mtrans() { \
02713                                                                          \
02714     shader_template<CLIP_VIEW_##mclipview,                               \
02715                     HEADLIGHT_##mhlight,                                 \
02716                     FOG_##mfog,                                          \
02717                     SHADOWS_##mshad,                                     \
02718                     AO_##mao,                                            \
02719                     OUTLINE_##moutl,                                     \
02720                     REFL_##mrefl,                                        \
02721                     TRANS_##mtrans >                                     \
02722                     (get_intersection_color<0>(),                        \
02723                      get_intersection_normal<0>(),                       \
02724                      Ka, Kd, Ks, phong_exp, Krefl,                       \
02725                      opacity, outline, outlinewidth, transmode);         \
02726   }
02727 
02728 
02729 //
02730 // Generate all of the 2^8 parameter combinations as 
02731 // template-specialized shaders...
02732 //
02733 
02734 DEFINE_CLOSEST_HIT(  on,  on,  on,  on,  on,  on,  on,  on )
02735 DEFINE_CLOSEST_HIT(  on,  on,  on,  on,  on,  on,  on, off )
02736 
02737 DEFINE_CLOSEST_HIT(  on,  on,  on,  on,  on,  on, off,  on )
02738 DEFINE_CLOSEST_HIT(  on,  on,  on,  on,  on,  on, off, off )
02739 
02740 DEFINE_CLOSEST_HIT(  on,  on,  on,  on,  on, off,  on,  on )
02741 DEFINE_CLOSEST_HIT(  on,  on,  on,  on,  on, off,  on, off )
02742 
02743 DEFINE_CLOSEST_HIT(  on,  on,  on,  on,  on, off, off,  on )
02744 DEFINE_CLOSEST_HIT(  on,  on,  on,  on,  on, off, off, off )
02745 
02746 DEFINE_CLOSEST_HIT(  on,  on,  on,  on, off,  on,  on,  on )
02747 DEFINE_CLOSEST_HIT(  on,  on,  on,  on, off,  on,  on, off )
02748 
02749 DEFINE_CLOSEST_HIT(  on,  on,  on,  on, off,  on, off,  on )
02750 DEFINE_CLOSEST_HIT(  on,  on,  on,  on, off,  on, off, off )
02751 
02752 DEFINE_CLOSEST_HIT(  on,  on,  on,  on, off, off,  on,  on )
02753 DEFINE_CLOSEST_HIT(  on,  on,  on,  on, off, off,  on, off )
02754 
02755 DEFINE_CLOSEST_HIT(  on,  on,  on,  on, off, off, off,  on )
02756 DEFINE_CLOSEST_HIT(  on,  on,  on,  on, off, off, off, off )
02757 
02758 
02759 DEFINE_CLOSEST_HIT(  on,  on,  on, off,  on,  on,  on,  on )
02760 DEFINE_CLOSEST_HIT(  on,  on,  on, off,  on,  on,  on, off )
02761 
02762 DEFINE_CLOSEST_HIT(  on,  on,  on, off,  on,  on, off,  on )
02763 DEFINE_CLOSEST_HIT(  on,  on,  on, off,  on,  on, off, off )
02764 
02765 DEFINE_CLOSEST_HIT(  on,  on,  on, off,  on, off,  on,  on )
02766 DEFINE_CLOSEST_HIT(  on,  on,  on, off,  on, off,  on, off )
02767 
02768 DEFINE_CLOSEST_HIT(  on,  on,  on, off,  on, off, off,  on )
02769 DEFINE_CLOSEST_HIT(  on,  on,  on, off,  on, off, off, off )
02770 
02771 DEFINE_CLOSEST_HIT(  on,  on,  on, off, off,  on,  on,  on )
02772 DEFINE_CLOSEST_HIT(  on,  on,  on, off, off,  on,  on, off )
02773 
02774 DEFINE_CLOSEST_HIT(  on,  on,  on, off, off,  on, off,  on )
02775 DEFINE_CLOSEST_HIT(  on,  on,  on, off, off,  on, off, off )
02776 
02777 DEFINE_CLOSEST_HIT(  on,  on,  on, off, off, off,  on,  on )
02778 DEFINE_CLOSEST_HIT(  on,  on,  on, off, off, off,  on, off )
02779 
02780 DEFINE_CLOSEST_HIT(  on,  on,  on, off, off, off, off,  on )
02781 DEFINE_CLOSEST_HIT(  on,  on,  on, off, off, off, off, off )
02782 
02783 
02784 
02785 DEFINE_CLOSEST_HIT(  on,  on, off,  on,  on,  on,  on,  on )
02786 DEFINE_CLOSEST_HIT(  on,  on, off,  on,  on,  on,  on, off )
02787 
02788 DEFINE_CLOSEST_HIT(  on,  on, off,  on,  on,  on, off,  on )
02789 DEFINE_CLOSEST_HIT(  on,  on, off,  on,  on,  on, off, off )
02790 
02791 DEFINE_CLOSEST_HIT(  on,  on, off,  on,  on, off,  on,  on )
02792 DEFINE_CLOSEST_HIT(  on,  on, off,  on,  on, off,  on, off )
02793 
02794 DEFINE_CLOSEST_HIT(  on,  on, off,  on,  on, off, off,  on )
02795 DEFINE_CLOSEST_HIT(  on,  on, off,  on,  on, off, off, off )
02796 
02797 DEFINE_CLOSEST_HIT(  on,  on, off,  on, off,  on,  on,  on )
02798 DEFINE_CLOSEST_HIT(  on,  on, off,  on, off,  on,  on, off )
02799 
02800 DEFINE_CLOSEST_HIT(  on,  on, off,  on, off,  on, off,  on )
02801 DEFINE_CLOSEST_HIT(  on,  on, off,  on, off,  on, off, off )
02802 
02803 DEFINE_CLOSEST_HIT(  on,  on, off,  on, off, off,  on,  on )
02804 DEFINE_CLOSEST_HIT(  on,  on, off,  on, off, off,  on, off )
02805 
02806 DEFINE_CLOSEST_HIT(  on,  on, off,  on, off, off, off,  on )
02807 DEFINE_CLOSEST_HIT(  on,  on, off,  on, off, off, off, off )
02808 
02809 
02810 DEFINE_CLOSEST_HIT(  on,  on, off, off,  on,  on,  on,  on )
02811 DEFINE_CLOSEST_HIT(  on,  on, off, off,  on,  on,  on, off )
02812 
02813 DEFINE_CLOSEST_HIT(  on,  on, off, off,  on,  on, off,  on )
02814 DEFINE_CLOSEST_HIT(  on,  on, off, off,  on,  on, off, off )
02815 
02816 DEFINE_CLOSEST_HIT(  on,  on, off, off,  on, off,  on,  on )
02817 DEFINE_CLOSEST_HIT(  on,  on, off, off,  on, off,  on, off )
02818 
02819 DEFINE_CLOSEST_HIT(  on,  on, off, off,  on, off, off,  on )
02820 DEFINE_CLOSEST_HIT(  on,  on, off, off,  on, off, off, off )
02821 
02822 DEFINE_CLOSEST_HIT(  on,  on, off, off, off,  on,  on,  on )
02823 DEFINE_CLOSEST_HIT(  on,  on, off, off, off,  on,  on, off )
02824 
02825 DEFINE_CLOSEST_HIT(  on,  on, off, off, off,  on, off,  on )
02826 DEFINE_CLOSEST_HIT(  on,  on, off, off, off,  on, off, off )
02827 
02828 DEFINE_CLOSEST_HIT(  on,  on, off, off, off, off,  on,  on )
02829 DEFINE_CLOSEST_HIT(  on,  on, off, off, off, off,  on, off )
02830 
02831 DEFINE_CLOSEST_HIT(  on,  on, off, off, off, off, off,  on )
02832 DEFINE_CLOSEST_HIT(  on,  on, off, off, off, off, off, off )
02833 
02834 //
02835 // block of 64 pgms
02836 //
02837 
02838 DEFINE_CLOSEST_HIT(  on, off,  on,  on,  on,  on,  on,  on )
02839 DEFINE_CLOSEST_HIT(  on, off,  on,  on,  on,  on,  on, off )
02840 
02841 DEFINE_CLOSEST_HIT(  on, off,  on,  on,  on,  on, off,  on )
02842 DEFINE_CLOSEST_HIT(  on, off,  on,  on,  on,  on, off, off )
02843 
02844 DEFINE_CLOSEST_HIT(  on, off,  on,  on,  on, off,  on,  on )
02845 DEFINE_CLOSEST_HIT(  on, off,  on,  on,  on, off,  on, off )
02846 
02847 DEFINE_CLOSEST_HIT(  on, off,  on,  on,  on, off, off,  on )
02848 DEFINE_CLOSEST_HIT(  on, off,  on,  on,  on, off, off, off )
02849 
02850 DEFINE_CLOSEST_HIT(  on, off,  on,  on, off,  on,  on,  on )
02851 DEFINE_CLOSEST_HIT(  on, off,  on,  on, off,  on,  on, off )
02852 
02853 DEFINE_CLOSEST_HIT(  on, off,  on,  on, off,  on, off,  on )
02854 DEFINE_CLOSEST_HIT(  on, off,  on,  on, off,  on, off, off )
02855 
02856 DEFINE_CLOSEST_HIT(  on, off,  on,  on, off, off,  on,  on )
02857 DEFINE_CLOSEST_HIT(  on, off,  on,  on, off, off,  on, off )
02858 
02859 DEFINE_CLOSEST_HIT(  on, off,  on,  on, off, off, off,  on )
02860 DEFINE_CLOSEST_HIT(  on, off,  on,  on, off, off, off, off )
02861 
02862 
02863 DEFINE_CLOSEST_HIT(  on, off,  on, off,  on,  on,  on,  on )
02864 DEFINE_CLOSEST_HIT(  on, off,  on, off,  on,  on,  on, off )
02865 
02866 DEFINE_CLOSEST_HIT(  on, off,  on, off,  on,  on, off,  on )
02867 DEFINE_CLOSEST_HIT(  on, off,  on, off,  on,  on, off, off )
02868 
02869 DEFINE_CLOSEST_HIT(  on, off,  on, off,  on, off,  on,  on )
02870 DEFINE_CLOSEST_HIT(  on, off,  on, off,  on, off,  on, off )
02871 
02872 DEFINE_CLOSEST_HIT(  on, off,  on, off,  on, off, off,  on )
02873 DEFINE_CLOSEST_HIT(  on, off,  on, off,  on, off, off, off )
02874 
02875 DEFINE_CLOSEST_HIT(  on, off,  on, off, off,  on,  on,  on )
02876 DEFINE_CLOSEST_HIT(  on, off,  on, off, off,  on,  on, off )
02877 
02878 DEFINE_CLOSEST_HIT(  on, off,  on, off, off,  on, off,  on )
02879 DEFINE_CLOSEST_HIT(  on, off,  on, off, off,  on, off, off )
02880 
02881 DEFINE_CLOSEST_HIT(  on, off,  on, off, off, off,  on,  on )
02882 DEFINE_CLOSEST_HIT(  on, off,  on, off, off, off,  on, off )
02883 
02884 DEFINE_CLOSEST_HIT(  on, off,  on, off, off, off, off,  on )
02885 DEFINE_CLOSEST_HIT(  on, off,  on, off, off, off, off, off )
02886 
02887 
02888 
02889 DEFINE_CLOSEST_HIT(  on, off, off,  on,  on,  on,  on,  on )
02890 DEFINE_CLOSEST_HIT(  on, off, off,  on,  on,  on,  on, off )
02891 
02892 DEFINE_CLOSEST_HIT(  on, off, off,  on,  on,  on, off,  on )
02893 DEFINE_CLOSEST_HIT(  on, off, off,  on,  on,  on, off, off )
02894 
02895 DEFINE_CLOSEST_HIT(  on, off, off,  on,  on, off,  on,  on )
02896 DEFINE_CLOSEST_HIT(  on, off, off,  on,  on, off,  on, off )
02897 
02898 DEFINE_CLOSEST_HIT(  on, off, off,  on,  on, off, off,  on )
02899 DEFINE_CLOSEST_HIT(  on, off, off,  on,  on, off, off, off )
02900 
02901 DEFINE_CLOSEST_HIT(  on, off, off,  on, off,  on,  on,  on )
02902 DEFINE_CLOSEST_HIT(  on, off, off,  on, off,  on,  on, off )
02903 
02904 DEFINE_CLOSEST_HIT(  on, off, off,  on, off,  on, off,  on )
02905 DEFINE_CLOSEST_HIT(  on, off, off,  on, off,  on, off, off )
02906 
02907 DEFINE_CLOSEST_HIT(  on, off, off,  on, off, off,  on,  on )
02908 DEFINE_CLOSEST_HIT(  on, off, off,  on, off, off,  on, off )
02909 
02910 DEFINE_CLOSEST_HIT(  on, off, off,  on, off, off, off,  on )
02911 DEFINE_CLOSEST_HIT(  on, off, off,  on, off, off, off, off )
02912 
02913 
02914 DEFINE_CLOSEST_HIT(  on, off, off, off,  on,  on,  on,  on )
02915 DEFINE_CLOSEST_HIT(  on, off, off, off,  on,  on,  on, off )
02916 
02917 DEFINE_CLOSEST_HIT(  on, off, off, off,  on,  on, off,  on )
02918 DEFINE_CLOSEST_HIT(  on, off, off, off,  on,  on, off, off )
02919 
02920 DEFINE_CLOSEST_HIT(  on, off, off, off,  on, off,  on,  on )
02921 DEFINE_CLOSEST_HIT(  on, off, off, off,  on, off,  on, off )
02922 
02923 DEFINE_CLOSEST_HIT(  on, off, off, off,  on, off, off,  on )
02924 DEFINE_CLOSEST_HIT(  on, off, off, off,  on, off, off, off )
02925 
02926 DEFINE_CLOSEST_HIT(  on, off, off, off, off,  on,  on,  on )
02927 DEFINE_CLOSEST_HIT(  on, off, off, off, off,  on,  on, off )
02928 
02929 DEFINE_CLOSEST_HIT(  on, off, off, off, off,  on, off,  on )
02930 DEFINE_CLOSEST_HIT(  on, off, off, off, off,  on, off, off )
02931 
02932 DEFINE_CLOSEST_HIT(  on, off, off, off, off, off,  on,  on )
02933 DEFINE_CLOSEST_HIT(  on, off, off, off, off, off,  on, off )
02934 
02935 DEFINE_CLOSEST_HIT(  on, off, off, off, off, off, off,  on )
02936 DEFINE_CLOSEST_HIT(  on, off, off, off, off, off, off, off )
02937 
02941 
02942 DEFINE_CLOSEST_HIT( off,  on,  on,  on,  on,  on,  on,  on )
02943 DEFINE_CLOSEST_HIT( off,  on,  on,  on,  on,  on,  on, off )
02944 
02945 DEFINE_CLOSEST_HIT( off,  on,  on,  on,  on,  on, off,  on )
02946 DEFINE_CLOSEST_HIT( off,  on,  on,  on,  on,  on, off, off )
02947 
02948 DEFINE_CLOSEST_HIT( off,  on,  on,  on,  on, off,  on,  on )
02949 DEFINE_CLOSEST_HIT( off,  on,  on,  on,  on, off,  on, off )
02950 
02951 DEFINE_CLOSEST_HIT( off,  on,  on,  on,  on, off, off,  on )
02952 DEFINE_CLOSEST_HIT( off,  on,  on,  on,  on, off, off, off )
02953 
02954 DEFINE_CLOSEST_HIT( off,  on,  on,  on, off,  on,  on,  on )
02955 DEFINE_CLOSEST_HIT( off,  on,  on,  on, off,  on,  on, off )
02956 
02957 DEFINE_CLOSEST_HIT( off,  on,  on,  on, off,  on, off,  on )
02958 DEFINE_CLOSEST_HIT( off,  on,  on,  on, off,  on, off, off )
02959 
02960 DEFINE_CLOSEST_HIT( off,  on,  on,  on, off, off,  on,  on )
02961 DEFINE_CLOSEST_HIT( off,  on,  on,  on, off, off,  on, off )
02962 
02963 DEFINE_CLOSEST_HIT( off,  on,  on,  on, off, off, off,  on )
02964 DEFINE_CLOSEST_HIT( off,  on,  on,  on, off, off, off, off )
02965 
02966 
02967 DEFINE_CLOSEST_HIT( off,  on,  on, off,  on,  on,  on,  on )
02968 DEFINE_CLOSEST_HIT( off,  on,  on, off,  on,  on,  on, off )
02969 
02970 DEFINE_CLOSEST_HIT( off,  on,  on, off,  on,  on, off,  on )
02971 DEFINE_CLOSEST_HIT( off,  on,  on, off,  on,  on, off, off )
02972 
02973 DEFINE_CLOSEST_HIT( off,  on,  on, off,  on, off,  on,  on )
02974 DEFINE_CLOSEST_HIT( off,  on,  on, off,  on, off,  on, off )
02975 
02976 DEFINE_CLOSEST_HIT( off,  on,  on, off,  on, off, off,  on )
02977 DEFINE_CLOSEST_HIT( off,  on,  on, off,  on, off, off, off )
02978 
02979 DEFINE_CLOSEST_HIT( off,  on,  on, off, off,  on,  on,  on )
02980 DEFINE_CLOSEST_HIT( off,  on,  on, off, off,  on,  on, off )
02981 
02982 DEFINE_CLOSEST_HIT( off,  on,  on, off, off,  on, off,  on )
02983 DEFINE_CLOSEST_HIT( off,  on,  on, off, off,  on, off, off )
02984 
02985 DEFINE_CLOSEST_HIT( off,  on,  on, off, off, off,  on,  on )
02986 DEFINE_CLOSEST_HIT( off,  on,  on, off, off, off,  on, off )
02987 
02988 DEFINE_CLOSEST_HIT( off,  on,  on, off, off, off, off,  on )
02989 DEFINE_CLOSEST_HIT( off,  on,  on, off, off, off, off, off )
02990 
02991 
02992 
02993 DEFINE_CLOSEST_HIT( off,  on, off,  on,  on,  on,  on,  on )
02994 DEFINE_CLOSEST_HIT( off,  on, off,  on,  on,  on,  on, off )
02995 
02996 DEFINE_CLOSEST_HIT( off,  on, off,  on,  on,  on, off,  on )
02997 DEFINE_CLOSEST_HIT( off,  on, off,  on,  on,  on, off, off )
02998 
02999 DEFINE_CLOSEST_HIT( off,  on, off,  on,  on, off,  on,  on )
03000 DEFINE_CLOSEST_HIT( off,  on, off,  on,  on, off,  on, off )
03001 
03002 DEFINE_CLOSEST_HIT( off,  on, off,  on,  on, off, off,  on )
03003 DEFINE_CLOSEST_HIT( off,  on, off,  on,  on, off, off, off )
03004 
03005 DEFINE_CLOSEST_HIT( off,  on, off,  on, off,  on,  on,  on )
03006 DEFINE_CLOSEST_HIT( off,  on, off,  on, off,  on,  on, off )
03007 
03008 DEFINE_CLOSEST_HIT( off,  on, off,  on, off,  on, off,  on )
03009 DEFINE_CLOSEST_HIT( off,  on, off,  on, off,  on, off, off )
03010 
03011 DEFINE_CLOSEST_HIT( off,  on, off,  on, off, off,  on,  on )
03012 DEFINE_CLOSEST_HIT( off,  on, off,  on, off, off,  on, off )
03013 
03014 DEFINE_CLOSEST_HIT( off,  on, off,  on, off, off, off,  on )
03015 DEFINE_CLOSEST_HIT( off,  on, off,  on, off, off, off, off )
03016 
03017 
03018 DEFINE_CLOSEST_HIT( off,  on, off, off,  on,  on,  on,  on )
03019 DEFINE_CLOSEST_HIT( off,  on, off, off,  on,  on,  on, off )
03020 
03021 DEFINE_CLOSEST_HIT( off,  on, off, off,  on,  on, off,  on )
03022 DEFINE_CLOSEST_HIT( off,  on, off, off,  on,  on, off, off )
03023 
03024 DEFINE_CLOSEST_HIT( off,  on, off, off,  on, off,  on,  on )
03025 DEFINE_CLOSEST_HIT( off,  on, off, off,  on, off,  on, off )
03026 
03027 DEFINE_CLOSEST_HIT( off,  on, off, off,  on, off, off,  on )
03028 DEFINE_CLOSEST_HIT( off,  on, off, off,  on, off, off, off )
03029 
03030 DEFINE_CLOSEST_HIT( off,  on, off, off, off,  on,  on,  on )
03031 DEFINE_CLOSEST_HIT( off,  on, off, off, off,  on,  on, off )
03032 
03033 DEFINE_CLOSEST_HIT( off,  on, off, off, off,  on, off,  on )
03034 DEFINE_CLOSEST_HIT( off,  on, off, off, off,  on, off, off )
03035 
03036 DEFINE_CLOSEST_HIT( off,  on, off, off, off, off,  on,  on )
03037 DEFINE_CLOSEST_HIT( off,  on, off, off, off, off,  on, off )
03038 
03039 DEFINE_CLOSEST_HIT( off,  on, off, off, off, off, off,  on )
03040 DEFINE_CLOSEST_HIT( off,  on, off, off, off, off, off, off )
03041 
03042 //
03043 // block of 64 pgms
03044 //
03045 
03046 DEFINE_CLOSEST_HIT( off, off,  on,  on,  on,  on,  on,  on )
03047 DEFINE_CLOSEST_HIT( off, off,  on,  on,  on,  on,  on, off )
03048 
03049 DEFINE_CLOSEST_HIT( off, off,  on,  on,  on,  on, off,  on )
03050 DEFINE_CLOSEST_HIT( off, off,  on,  on,  on,  on, off, off )
03051 
03052 DEFINE_CLOSEST_HIT( off, off,  on,  on,  on, off,  on,  on )
03053 DEFINE_CLOSEST_HIT( off, off,  on,  on,  on, off,  on, off )
03054 
03055 DEFINE_CLOSEST_HIT( off, off,  on,  on,  on, off, off,  on )
03056 DEFINE_CLOSEST_HIT( off, off,  on,  on,  on, off, off, off )
03057 
03058 DEFINE_CLOSEST_HIT( off, off,  on,  on, off,  on,  on,  on )
03059 DEFINE_CLOSEST_HIT( off, off,  on,  on, off,  on,  on, off )
03060 
03061 DEFINE_CLOSEST_HIT( off, off,  on,  on, off,  on, off,  on )
03062 DEFINE_CLOSEST_HIT( off, off,  on,  on, off,  on, off, off )
03063 
03064 DEFINE_CLOSEST_HIT( off, off,  on,  on, off, off,  on,  on )
03065 DEFINE_CLOSEST_HIT( off, off,  on,  on, off, off,  on, off )
03066 
03067 DEFINE_CLOSEST_HIT( off, off,  on,  on, off, off, off,  on )
03068 DEFINE_CLOSEST_HIT( off, off,  on,  on, off, off, off, off )
03069 
03070 
03071 DEFINE_CLOSEST_HIT( off, off,  on, off,  on,  on,  on,  on )
03072 DEFINE_CLOSEST_HIT( off, off,  on, off,  on,  on,  on, off )
03073 
03074 DEFINE_CLOSEST_HIT( off, off,  on, off,  on,  on, off,  on )
03075 DEFINE_CLOSEST_HIT( off, off,  on, off,  on,  on, off, off )
03076 
03077 DEFINE_CLOSEST_HIT( off, off,  on, off,  on, off,  on,  on )
03078 DEFINE_CLOSEST_HIT( off, off,  on, off,  on, off,  on, off )
03079 
03080 DEFINE_CLOSEST_HIT( off, off,  on, off,  on, off, off,  on )
03081 DEFINE_CLOSEST_HIT( off, off,  on, off,  on, off, off, off )
03082 
03083 DEFINE_CLOSEST_HIT( off, off,  on, off, off,  on,  on,  on )
03084 DEFINE_CLOSEST_HIT( off, off,  on, off, off,  on,  on, off )
03085 
03086 DEFINE_CLOSEST_HIT( off, off,  on, off, off,  on, off,  on )
03087 DEFINE_CLOSEST_HIT( off, off,  on, off, off,  on, off, off )
03088 
03089 DEFINE_CLOSEST_HIT( off, off,  on, off, off, off,  on,  on )
03090 DEFINE_CLOSEST_HIT( off, off,  on, off, off, off,  on, off )
03091 
03092 DEFINE_CLOSEST_HIT( off, off,  on, off, off, off, off,  on )
03093 DEFINE_CLOSEST_HIT( off, off,  on, off, off, off, off, off )
03094 
03095 
03096 
03097 DEFINE_CLOSEST_HIT( off, off, off,  on,  on,  on,  on,  on )
03098 DEFINE_CLOSEST_HIT( off, off, off,  on,  on,  on,  on, off )
03099 
03100 DEFINE_CLOSEST_HIT( off, off, off,  on,  on,  on, off,  on )
03101 DEFINE_CLOSEST_HIT( off, off, off,  on,  on,  on, off, off )
03102 
03103 DEFINE_CLOSEST_HIT( off, off, off,  on,  on, off,  on,  on )
03104 DEFINE_CLOSEST_HIT( off, off, off,  on,  on, off,  on, off )
03105 
03106 DEFINE_CLOSEST_HIT( off, off, off,  on,  on, off, off,  on )
03107 DEFINE_CLOSEST_HIT( off, off, off,  on,  on, off, off, off )
03108 
03109 DEFINE_CLOSEST_HIT( off, off, off,  on, off,  on,  on,  on )
03110 DEFINE_CLOSEST_HIT( off, off, off,  on, off,  on,  on, off )
03111 
03112 DEFINE_CLOSEST_HIT( off, off, off,  on, off,  on, off,  on )
03113 DEFINE_CLOSEST_HIT( off, off, off,  on, off,  on, off, off )
03114 
03115 DEFINE_CLOSEST_HIT( off, off, off,  on, off, off,  on,  on )
03116 DEFINE_CLOSEST_HIT( off, off, off,  on, off, off,  on, off )
03117 
03118 DEFINE_CLOSEST_HIT( off, off, off,  on, off, off, off,  on )
03119 DEFINE_CLOSEST_HIT( off, off, off,  on, off, off, off, off )
03120 
03121 
03122 DEFINE_CLOSEST_HIT( off, off, off, off,  on,  on,  on,  on )
03123 DEFINE_CLOSEST_HIT( off, off, off, off,  on,  on,  on, off )
03124 
03125 DEFINE_CLOSEST_HIT( off, off, off, off,  on,  on, off,  on )
03126 DEFINE_CLOSEST_HIT( off, off, off, off,  on,  on, off, off )
03127 
03128 DEFINE_CLOSEST_HIT( off, off, off, off,  on, off,  on,  on )
03129 DEFINE_CLOSEST_HIT( off, off, off, off,  on, off,  on, off )
03130 
03131 DEFINE_CLOSEST_HIT( off, off, off, off,  on, off, off,  on )
03132 DEFINE_CLOSEST_HIT( off, off, off, off,  on, off, off, off )
03133 
03134 DEFINE_CLOSEST_HIT( off, off, off, off, off,  on,  on,  on )
03135 DEFINE_CLOSEST_HIT( off, off, off, off, off,  on,  on, off )
03136 
03137 DEFINE_CLOSEST_HIT( off, off, off, off, off,  on, off,  on )
03138 DEFINE_CLOSEST_HIT( off, off, off, off, off,  on, off, off )
03139 
03140 DEFINE_CLOSEST_HIT( off, off, off, off, off, off,  on,  on )
03141 DEFINE_CLOSEST_HIT( off, off, off, off, off, off,  on, off )
03142 
03143 DEFINE_CLOSEST_HIT( off, off, off, off, off, off, off,  on )
03144 DEFINE_CLOSEST_HIT( off, off, off, off, off, off, off, off )
03145 
03146 
03147 #undef CLIP_VIEW_on
03148 #undef CLIP_VIEW_off
03149 #undef HEADLIGHT_on
03150 #undef HEADLIGHT_off
03151 #undef FOG_on
03152 #undef FOG_off
03153 #undef SHADOWS_on
03154 #undef SHADOWS_off
03155 #undef AO_on
03156 #undef AO_off
03157 #undef OUTLINE_on
03158 #undef OUTLINE_off
03159 #undef REFL_on       
03160 #undef REFL_off           
03161 #undef TRANS_on         
03162 #undef TRANS_off          
03163 #undef DEFINE_CLOSEST_HIT
03164 
03165 #endif // ORT_USE_TEMPLATE_SHADERS
03166 
03167 

Generated on Tue Oct 15 02:44:50 2024 for VMD (current) by doxygen1.2.14 written by Dimitri van Heesch, © 1997-2002