1.3.49
 
Loading...
Searching...
No Matches
rayHit.cu
Go to the documentation of this file.
1
17#include <optix_world.h>
18#include <optixu/optixu_math_namespace.h>
19#include <optixu/optixu_matrix_namespace.h>
20
21#include "RayTracing.cuh"
22
23using namespace optix;
24
25rtDeclareVariable(optix::Ray, ray, rtCurrentRay, );
26rtDeclareVariable(float, t_hit, rtIntersectionDistance, );
27rtDeclareVariable(PerRayData, prd, rtPayload, );
28
29rtDeclareVariable(unsigned int, UUID, attribute UUID, );
30
31RT_PROGRAM void closest_hit_direct() {
32
33 uint objID = objectID[UUID];
34
35 if ((periodic_flag.x == 1 || periodic_flag.y == 1) && primitive_type[objID] == 5) { // periodic boundary condition
36
37 prd.hit_periodic_boundary = true;
38
39 float3 ray_origin = ray.origin + t_hit * ray.direction;
40
41 float eps = 1e-5;
42
43 float2 xbounds = make_float2(bbox_vertices[make_uint2(0, 0)].x, bbox_vertices[make_uint2(1, 1)].x);
44 float2 ybounds = make_float2(bbox_vertices[make_uint2(0, 0)].y, bbox_vertices[make_uint2(1, 1)].y);
45
46 float width_x = xbounds.y - xbounds.x;
47 float width_y = ybounds.y - ybounds.x;
48
49 prd.periodic_hit = ray_origin;
50 if (periodic_flag.x == 1 && fabs(ray_origin.x - xbounds.x) <= eps) { //-x facing boundary
51 prd.periodic_hit.x += +width_x - eps;
52 } else if (periodic_flag.x == 1 && fabs(ray_origin.x - xbounds.y) <= eps) { //+x facing boundary
53 prd.periodic_hit.x += -width_x + eps;
54 } else if (periodic_flag.y == 1 && fabs(ray_origin.y - ybounds.x) <= eps) { //-y facing boundary
55 prd.periodic_hit.y += +width_y - eps;
56 } else if (periodic_flag.y == 1 && fabs(ray_origin.y - ybounds.y) <= eps) { //+y facing boundary
57 prd.periodic_hit.y += -width_y + eps;
58 }
59 }
60};
61
62RT_PROGRAM void closest_hit_diffuse() {
63
64 uint origin_UUID = prd.origin_UUID;
65
66 uint objID = objectID[UUID];
67
68 if ((periodic_flag.x == 1 || periodic_flag.y == 1) && primitive_type[objID] == 5) { // periodic boundary condition
69
70 prd.hit_periodic_boundary = true;
71
72 float3 ray_origin = ray.origin + t_hit * ray.direction;
73
74 float eps = 1e-5;
75
76 float2 xbounds = make_float2(bbox_vertices[make_uint2(0, 0)].x, bbox_vertices[make_uint2(1, 1)].x);
77 float2 ybounds = make_float2(bbox_vertices[make_uint2(0, 0)].y, bbox_vertices[make_uint2(1, 1)].y);
78
79 float width_x = xbounds.y - xbounds.x;
80 float width_y = ybounds.y - ybounds.x;
81
82 prd.periodic_hit = ray_origin;
83 if (periodic_flag.x == 1 && fabs(ray_origin.x - xbounds.x) <= eps) { //-x facing boundary
84 prd.periodic_hit.x += +width_x - eps;
85 } else if (periodic_flag.x == 1 && fabs(ray_origin.x - xbounds.y) <= eps) { //+x facing boundary
86 prd.periodic_hit.x += -width_x + eps;
87 } else if (periodic_flag.y == 1 && fabs(ray_origin.y - ybounds.x) <= eps) { //-y facing boundary
88 prd.periodic_hit.y += +width_y - eps;
89 } else if (periodic_flag.y == 1 && fabs(ray_origin.y - ybounds.y) <= eps) { //+y facing boundary
90 prd.periodic_hit.y += -width_y + eps;
91 }
92
93 } else {
94
95 // Note: UUID corresponds to the object that the ray hit (i.e., where energy is coming from), and UUID_origin is the object the ray originated from (i.e., where the energy is being recieved)
96
97 // find out if we hit top or bottom surface
98 float3 normal;
99
100 float m[16];
101 for (uint i = 0; i < 16; i++) {
102 m[i] = transform_matrix[optix::make_uint2(i, objID)];
103 }
104
105 if (primitive_type[objID] == 0 || primitive_type[objID] == 3) { // hit patch or tile
106 float3 s0 = make_float3(0, 0, 0);
107 float3 s1 = make_float3(1, 0, 0);
108 float3 s2 = make_float3(0, 1, 0);
109 d_transformPoint(m, s0);
110 d_transformPoint(m, s1);
111 d_transformPoint(m, s2);
112 normal = cross(s1 - s0, s2 - s0);
113 } else if (primitive_type[UUID] == 1) { // hit triangle
114 float3 v0 = make_float3(0, 0, 0);
115 d_transformPoint(m, v0);
116 float3 v1 = make_float3(0, 1, 0);
117 d_transformPoint(m, v1);
118 float3 v2 = make_float3(1, 1, 0);
119 d_transformPoint(m, v2);
120 normal = cross(v1 - v0, v2 - v1);
121 } else if (primitive_type[UUID] == 2) { // hit disk
122 float3 v0 = make_float3(0, 0, 0);
123 d_transformPoint(m, v0);
124 float3 v1 = make_float3(1, 0, 0);
125 d_transformPoint(m, v1);
126 float3 v2 = make_float3(0, 1, 0);
127 d_transformPoint(m, v2);
128 normal = cross(v1 - v0, v2 - v0);
129 } else if (primitive_type[UUID] == 4) { // hit voxel
130 float3 vmin = make_float3(-0.5, -0.5, -0.5);
131 d_transformPoint(m, vmin);
132 float3 vmax = make_float3(0.5, 0.5, 0.5);
133 d_transformPoint(m, vmax);
134 }
135 normal = normalize(normal);
136
137 bool face = dot(normal, ray.direction) < 0;
138
139 int b = -1;
140 for (int b_global = 0; b_global < Nbands_global; b_global++) {
141
142 if (band_launch_flag[b_global] == 0) {
143 continue;
144 }
145 b++;
146
147 size_t ind_origin = Nbands_launch * origin_UUID + b;
148 size_t ind_hit = Nbands_launch * UUID + b;
149
150 double strength;
151 if (face || primitive_type[objID] == 4) {
152 strength = radiation_out_top[ind_hit] * prd.strength;
153 } else {
154 strength = radiation_out_bottom[ind_hit] * prd.strength;
155 }
156
157 if (strength == 0) {
158 continue;
159 }
160
161 size_t radprop_ind_global = Nprimitives * Nbands_global * prd.source_ID + Nbands_global * origin_UUID + b_global;
162 float t_rho = rho[radprop_ind_global];
163 float t_tau = tau[radprop_ind_global];
164
165 if (primitive_type[objectID[origin_UUID]] == 4) { // ray was launched from voxel
166
167 // float kappa = t_rho; //just a reminder that rho is actually the absorption coefficient
168 // float sigma_s = t_tau; //just a reminder that tau is actually the scattering coefficient
169 // float beta = kappa + sigma_s;
170 //
171 // // absorption
172 // atomicAdd(&radiation_in[ind_origin], strength * exp(-beta * 0.5 * prd.area) * kappa / beta);
173 //
174 // // scattering
175 // atomicAdd(&scatter_buff_top[ind_origin], strength * exp(-beta * 0.5 * prd.area) * sigma_s / beta);
176
177 } else { // ray was NOT launched from voxel
178
179 // absorption
180 atomicAdd(&radiation_in[ind_origin], strength * (1.f - t_rho - t_tau));
181
182 if ((t_rho > 0 || t_tau > 0) && strength > 0) {
183 if (prd.face) { // reflection from top, transmission from bottom
184 atomicFloatAdd(&scatter_buff_top[ind_origin], strength * t_rho); // reflection
185 atomicFloatAdd(&scatter_buff_bottom[ind_origin], strength * t_tau); // transmission
186 } else { // reflection from bottom, transmission from top
187 atomicFloatAdd(&scatter_buff_bottom[ind_origin], strength * t_rho); // reflection
188 atomicFloatAdd(&scatter_buff_top[ind_origin], strength * t_tau); // transmission
189 }
190 }
191 if (Ncameras > 0) {
192 size_t indc = prd.source_ID * Nprimitives * Nbands_global * Ncameras + origin_UUID * Nbands_global * Ncameras + b_global * Ncameras + camera_ID;
193 float t_rho_cam = rho_cam[indc];
194 float t_tau_cam = tau_cam[indc];
195 if ((t_rho_cam > 0 || t_tau_cam > 0) && strength > 0) {
196 if (prd.face) { // reflection from top, transmission from bottom
197 atomicFloatAdd(&scatter_buff_top_cam[ind_origin], strength * t_rho_cam); // reflection
198 atomicFloatAdd(&scatter_buff_bottom_cam[ind_origin], strength * t_tau_cam); // transmission
199 } else { // reflection from bottom, transmission from top
200 atomicFloatAdd(&scatter_buff_bottom_cam[ind_origin], strength * t_rho_cam); // reflection
201 atomicFloatAdd(&scatter_buff_top_cam[ind_origin], strength * t_tau_cam); // transmission
202 }
203 }
204 }
205 }
206
207 // if( primitive_type[UUID] == 4 ){ //if we hit a voxel, reduce strength and launch another ray
208 // optix::Ray ray_transmit = optix::make_Ray(ray.origin+(t_hit+prd.area+1e-5)*ray.direction, ray.direction, ray.ray_type, 1e-4, RT_DEFAULT_MAX);
209 // PerRayData prd_transmit = prd;
210 // float beta = rho[UUID]+tau[UUID];
211 // prd_transmit.strength = prd.strength*(1.f-exp(-beta*0.5*prd.area));
212 // rtTrace( top_object, ray_transmit, prd_transmit);
213 // }
214 }
215 }
216}
217
218RT_PROGRAM void closest_hit_camera() {
219
220 uint objID = objectID[UUID];
221
222 if ((periodic_flag.x == 1 || periodic_flag.y == 1) && primitive_type[objID] == 5) { // periodic boundary condition
223
224 prd.hit_periodic_boundary = true;
225
226 float3 ray_origin = ray.origin + t_hit * ray.direction;
227
228 float eps = 1e-5;
229
230 float2 xbounds = make_float2(bbox_vertices[make_uint2(0, 0)].x, bbox_vertices[make_uint2(1, 1)].x);
231 float2 ybounds = make_float2(bbox_vertices[make_uint2(0, 0)].y, bbox_vertices[make_uint2(1, 1)].y);
232
233 float width_x = xbounds.y - xbounds.x;
234 float width_y = ybounds.y - ybounds.x;
235
236 prd.periodic_hit = ray_origin;
237 if (periodic_flag.x == 1 && fabs(ray_origin.x - xbounds.x) <= eps) { //-x facing boundary
238 prd.periodic_hit.x += +width_x - eps;
239 } else if (periodic_flag.x == 1 && fabs(ray_origin.x - xbounds.y) <= eps) { //+x facing boundary
240 prd.periodic_hit.x += -width_x + eps;
241 } else if (periodic_flag.y == 1 && fabs(ray_origin.y - ybounds.x) <= eps) { //-y facing boundary
242 prd.periodic_hit.y += +width_y - eps;
243 } else if (periodic_flag.y == 1 && fabs(ray_origin.y - ybounds.y) <= eps) { //+y facing boundary
244 prd.periodic_hit.y += -width_y + eps;
245 }
246
247 } else {
248
249 // Note: UUID corresponds to the object that the ray hit (i.e., where energy is coming from), and UUID_origin is the object the ray originated from (i.e., where the energy is being received)
250
251 // find out if we hit top or bottom surface per each ray
252 float3 normal;
253
254 float m[16];
255 for (uint i = 0; i < 16; i++) {
256 m[i] = transform_matrix[optix::make_uint2(i, objID)];
257 }
258
259 if (primitive_type[objID] == 0 || primitive_type[objID] == 3) { // hit patch or tile
260 float3 s0 = make_float3(0, 0, 0);
261 float3 s1 = make_float3(1, 0, 0);
262 float3 s2 = make_float3(0, 1, 0);
263 d_transformPoint(m, s0);
264 d_transformPoint(m, s1);
265 d_transformPoint(m, s2);
266 normal = cross(s1 - s0, s2 - s0);
267 } else if (primitive_type[UUID] == 1) { // hit triangle
268 float3 v0 = make_float3(0, 0, 0);
269 d_transformPoint(m, v0);
270 float3 v1 = make_float3(0, 1, 0);
271 d_transformPoint(m, v1);
272 float3 v2 = make_float3(1, 1, 0);
273 d_transformPoint(m, v2);
274 normal = cross(v1 - v0, v2 - v1);
275 } else if (primitive_type[UUID] == 2) { // hit disk
276 float3 v0 = make_float3(0, 0, 0);
277 d_transformPoint(m, v0);
278 float3 v1 = make_float3(1, 0, 0);
279 d_transformPoint(m, v1);
280 float3 v2 = make_float3(0, 1, 0);
281 d_transformPoint(m, v2);
282 normal = cross(v1 - v0, v2 - v0);
283 } else if (primitive_type[UUID] == 4) { // hit voxel
284 float3 vmin = make_float3(-0.5, -0.5, -0.5);
285 d_transformPoint(m, vmin);
286 float3 vmax = make_float3(0.5, 0.5, 0.5);
287 d_transformPoint(m, vmax);
288 }
289 normal = normalize(normal);
290
291 bool face = dot(normal, ray.direction) < 0;
292
293
294 float3 camera_normal = d_rotatePoint(make_float3(0, 0, 1), -0.5 * M_PI + camera_direction.x, 0.5f * M_PI - camera_direction.y);
295
296 double strength;
297 for (size_t b = 0; b < Nbands_launch; b++) {
298
299 if (face || primitive_type[objID] == 4) {
300 strength = radiation_out_top[Nbands_launch * UUID + b] * prd.strength; // this one /fabs(dot())
301 } else {
302 strength = radiation_out_bottom[Nbands_launch * UUID + b] * prd.strength;
303 }
304
305 // specular reflection
306
307 double strength_spec = 0;
308 if (specular_reflection_enabled > 0 && specular_exponent[objID] > 0.f) {
309 for (int rr = 0; rr < Nsources; rr++) {
310
311 // light direction
312 float3 light_direction;
313 float spec = 0;
314 if (source_types[rr] == 0 || source_types[rr] == 2) { // collimated source or sunsphere source
315
316 light_direction = normalize(source_positions[rr]);
317 if (face) {
318 spec = radiation_out_top[Nbands_launch * UUID + b];
319 } else {
320 spec = radiation_out_bottom[Nbands_launch * UUID + b];
321 }
322
323 } else { // sphere, disk or rectangular source
324 //\todo Need to add generic functions to sample points on sphere, disk and rectangle source surfaces
325 }
326
327 // float3 specular_direction = normalize(2 * fabs(dot(light_direction, normal)) * normal - light_direction);
328 float3 specular_direction = normalize(light_direction - ray.direction);
329
330 // specular reflection
331 float exponent = specular_exponent[objID];
332 double scale_coefficient = 1.0;
333 if (specular_reflection_enabled == 2) { // if we are using the scale coefficient
334 scale_coefficient = specular_scale[objID];
335 }
336 strength_spec += spec * scale_coefficient * pow(max(0.f, dot(specular_direction, normal)), exponent) * (exponent + 2.f) /
337 (double(launch_dim.x) * 2.f * M_PI); // launch_dim.x is the number of rays launched per pixel, so we divide by it to get the average flux per ray. (exponent+2)/2pi normalizes reflected distribution to unity.
338 }
339 }
340
341 // absorption
342
343 atomicAdd(&radiation_in_camera[Nbands_launch * prd.origin_UUID + b], strength + strength_spec);
344 }
345 }
346}
347
348RT_PROGRAM void closest_hit_pixel_label() {
349
350 uint origin_UUID = prd.origin_UUID;
351
352 uint objID = objectID[UUID];
353
354 if ((periodic_flag.x == 1 || periodic_flag.y == 1) && primitive_type[objID] == 5) { // periodic boundary condition
355
356 prd.hit_periodic_boundary = true;
357
358 float3 ray_origin = ray.origin + t_hit * ray.direction;
359
360 float eps = 1e-5;
361
362 float2 xbounds = make_float2(bbox_vertices[make_uint2(0, 0)].x, bbox_vertices[make_uint2(1, 1)].x);
363 float2 ybounds = make_float2(bbox_vertices[make_uint2(0, 0)].y, bbox_vertices[make_uint2(1, 1)].y);
364
365 float width_x = xbounds.y - xbounds.x;
366 float width_y = ybounds.y - ybounds.x;
367
368 prd.periodic_hit = ray_origin;
369 if (periodic_flag.x == 1 && fabs(ray_origin.x - xbounds.x) <= eps) { //-x facing boundary
370 prd.periodic_hit.x += +width_x - eps;
371 } else if (periodic_flag.x == 1 && fabs(ray_origin.x - xbounds.y) <= eps) { //+x facing boundary
372 prd.periodic_hit.x += -width_x + eps;
373 } else if (periodic_flag.y == 1 && fabs(ray_origin.y - ybounds.x) <= eps) { //-y facing boundary
374 prd.periodic_hit.y += +width_y - eps;
375 } else if (periodic_flag.y == 1 && fabs(ray_origin.y - ybounds.y) <= eps) { //+y facing boundary
376 prd.periodic_hit.y += -width_y + eps;
377 }
378
379 } else {
380
381 // Note: UUID corresponds to the object that the ray hit (i.e., where energy is coming from), and UUID_origin is the object the ray originated from (i.e., where the energy is being received)
382 // Note: We are reserving a value of 0 for the sky, so we will store UUID+1
383 camera_pixel_label[origin_UUID] = UUID + 1;
384
385 float depth = prd.strength + t_hit;
386 float3 camera_direction3 = d_rotatePoint(make_float3(1, 0, 0), -0.5 * M_PI + camera_direction.x, 0.5f * M_PI - camera_direction.y);
387 camera_pixel_depth[origin_UUID] = abs(dot(camera_direction3, ray.direction)) * depth;
388 }
389}
390
391RT_PROGRAM void miss_direct() {
392
393 uint objID = objectID[prd.origin_UUID];
394
395 int b = -1;
396 for (int b_global = 0; b_global < Nbands_global; b_global++) {
397
398 if (band_launch_flag[b_global] == 0) {
399 continue;
400 }
401 b++;
402
403 size_t ind_origin = Nbands_launch * prd.origin_UUID + b;
404
405 size_t radprop_ind_global = Nprimitives * Nbands_global * prd.source_ID + Nbands_global * prd.origin_UUID + b_global;
406 float t_rho = rho[radprop_ind_global];
407 float t_tau = tau[radprop_ind_global];
408
409 double strength = prd.strength * source_fluxes[prd.source_ID * Nbands_launch + b];
410
411 // absorption
412 atomicAdd(&radiation_in[ind_origin], strength * (1.f - t_rho - t_tau));
413
414 if (t_rho > 0 || t_tau > 0) {
415 if (prd.face) { // reflection from top, transmission from bottom
416 atomicFloatAdd(&scatter_buff_top[ind_origin], strength * t_rho); // reflection
417 atomicFloatAdd(&scatter_buff_bottom[ind_origin], strength * t_tau); // transmission
418 } else { // reflection from bottom, transmission from top
419 atomicFloatAdd(&scatter_buff_bottom[ind_origin], strength * t_rho); // reflection
420 atomicFloatAdd(&scatter_buff_top[ind_origin], strength * t_tau); // transmission
421 }
422 }
423 if (Ncameras > 0) {
424 size_t indc = prd.source_ID * Nprimitives * Nbands_global * Ncameras + prd.origin_UUID * Nbands_global * Ncameras + b_global * Ncameras + camera_ID;
425 float t_rho_cam = rho_cam[indc];
426 float t_tau_cam = tau_cam[indc];
427 if ((t_rho_cam > 0 || t_tau_cam > 0) && strength > 0) {
428 if (prd.face) { // reflection from top, transmission from bottom
429 atomicFloatAdd(&scatter_buff_top_cam[ind_origin], strength * t_rho_cam); // reflection
430 atomicFloatAdd(&scatter_buff_bottom_cam[ind_origin], strength * t_tau_cam); // transmission
431 } else { // reflection from bottom, transmission from top
432 atomicFloatAdd(&scatter_buff_bottom_cam[ind_origin], strength * t_rho_cam); // reflection
433 atomicFloatAdd(&scatter_buff_top_cam[ind_origin], strength * t_tau_cam); // transmission
434 }
435 }
436 }
437 }
438}
439
440RT_PROGRAM void miss_diffuse() {
441
442 // double strength;
443 // if (prd.face || primitive_type[objectID[prd.origin_UUID]] == 3) {
444 // strength = radiation_out_top[prd.origin_UUID] * prd.strength * prd.area;
445 // } else {
446 // strength = radiation_out_bottom[prd.origin_UUID] * prd.strength * prd.area;
447 // }
448 //
449 // atomicFloatAdd(&Rsky[prd.origin_UUID], strength);
450
451 int b = -1;
452 for (size_t b_global = 0; b_global < Nbands_global; b_global++) {
453
454 if (band_launch_flag[b_global] == 0) {
455 continue;
456 }
457 b++;
458
459 if (diffuse_flux[b] > 0.f) {
460
461 size_t ind_origin = Nbands_launch * prd.origin_UUID + b;
462
463 size_t radprop_ind_global = Nprimitives * Nbands_global * prd.source_ID + Nbands_global * prd.origin_UUID + b_global;
464 float t_rho = rho[radprop_ind_global];
465 float t_tau = tau[radprop_ind_global];
466
467 if (primitive_type[objectID[prd.origin_UUID]] == 4) { // ray was launched from voxel
468
469 float kappa = t_rho; // just a reminder that rho is actually the absorption coefficient
470 float sigma_s = t_tau; // just a reminder that tau is actually the scattering coefficient
471 float beta = kappa + sigma_s;
472
473 // absorption
474 atomicAdd(&radiation_in[ind_origin], diffuse_flux[b] * prd.strength * kappa / beta);
475
476 // scattering
477 atomicAdd(&scatter_buff_top[ind_origin], diffuse_flux[b] * prd.strength * sigma_s / beta);
478
479 } else { // ray was NOT launched from voxel
480
481 float fd = 1.f;
482 if (diffuse_extinction[b] > 0.f) {
483 float psi = acos_safe(dot(diffuse_peak_dir[b], ray.direction));
484 if (psi < M_PI / 180.f) {
485 fd = powf(M_PI / 180.f, -diffuse_extinction[b]) * diffuse_dist_norm[b]; // Replace 'pow' by 'powf' in
486 } else {
487 fd = powf(psi, -diffuse_extinction[b]) * diffuse_dist_norm[b];
488 }
489 }
490
491 float strength = fd * diffuse_flux[b] * prd.strength;
492
493 // absorption
494 atomicAdd(&radiation_in[ind_origin], strength * (1.f - t_rho - t_tau));
495
496 if (t_rho > 0 || t_tau > 0) {
497 if (prd.face) { // reflection from top, transmission from bottom
498 atomicFloatAdd(&scatter_buff_top[ind_origin], strength * t_rho); // reflection
499 atomicFloatAdd(&scatter_buff_bottom[ind_origin], strength * t_tau); // transmission
500 } else { // reflection from bottom, transmission from top
501 atomicFloatAdd(&scatter_buff_bottom[ind_origin], strength * t_rho); // reflection
502 atomicFloatAdd(&scatter_buff_top[ind_origin], strength * t_tau); // transmission
503 }
504 }
505 if (Ncameras > 0) {
506 size_t indc = prd.source_ID * Nprimitives * Nbands_global * Ncameras + prd.origin_UUID * Nbands_global * Ncameras + b_global * Ncameras + camera_ID;
507 float t_rho_cam = rho_cam[indc];
508 float t_tau_cam = tau_cam[indc];
509 if ((t_rho_cam > 0 || t_tau_cam > 0) && prd.strength > 0) {
510 if (prd.face) { // reflection from top, transmission from bottom
511 atomicFloatAdd(&scatter_buff_top_cam[ind_origin], strength * t_rho_cam); // reflection
512 atomicFloatAdd(&scatter_buff_bottom_cam[ind_origin], strength * t_tau_cam); // transmission
513 } else { // reflection from bottom, transmission from top
514 atomicFloatAdd(&scatter_buff_bottom_cam[ind_origin], strength * t_rho_cam); // reflection
515 atomicFloatAdd(&scatter_buff_top_cam[ind_origin], strength * t_tau_cam); // transmission
516 }
517 }
518 }
519 }
520 }
521 }
522}
523
524RT_PROGRAM void miss_camera() {
525
526 for (size_t b = 0; b < Nbands_launch; b++) {
527
528 if (diffuse_flux[b] > 0.f) {
529
530 float fd = 1.f;
531 if (diffuse_extinction[b] > 0.f) {
532 float psi = acos_safe(dot(diffuse_peak_dir[b], ray.direction));
533 if (psi < M_PI / 180.f) {
534 fd = powf(float(M_PI) / 180.f, -diffuse_extinction[b]) * diffuse_dist_norm[b];
535 } else {
536 fd = powf(psi, -diffuse_extinction[b]) * diffuse_dist_norm[b];
537 }
538 }
539
540 // absorption
541 atomicAdd(&radiation_in_camera[Nbands_launch * prd.origin_UUID + b], fd * diffuse_flux[b] * prd.strength);
542 }
543 }
544}
545
546RT_PROGRAM void miss_pixel_label() {
547
548 camera_pixel_depth[prd.origin_UUID] = -1;
549}