1.3.49
 
Loading...
Searching...
No Matches
rayGeneration.cu
Go to the documentation of this file.
1
16#include <optix_world.h>
17#include "RayTracing.cuh"
18
19using namespace optix;
20
21RT_PROGRAM void direct_raygen() {
22
23 uint Nrays = launch_dim.x * launch_dim.y;
24 uint ray_index = launch_dim.x * launch_index.y + launch_index.x;
25
26 PerRayData prd;
27 prd.seed = tea<16>(ray_index + Nrays * launch_index.z, random_seed);
28
29 uint objID = launch_offset + launch_index.z;
30
31 uint pID = primitiveID[objID];
32 uint ptype = primitive_type[objID];
33 int puvID = uvID[objID];
34
35 float3 sp;
36
37 float3 normal;
38
39 // transformation matrix
40 float m_trans[16];
41 for (uint i = 0; i < 16; i++) {
42 m_trans[i] = transform_matrix[optix::make_uint2(i, objID)];
43 }
44
45 // looping over sub-patches
46 int NX = object_subdivisions[objID].x;
47 int NY = object_subdivisions[objID].y;
48 for (int jj = 0; jj < NY; jj++) {
49 for (int ii = 0; ii < NX; ii++) {
50
51 uint UUID = pID + jj * NX + ii;
52
53 // two random samples [0,1]
54 float Rx = rnd(prd.seed);
55 float Ry = rnd(prd.seed);
56
57 if (ptype == 0 || ptype == 3) { // Patch or Tile
58
59 uint Nx = launch_dim.x;
60 uint Ny = launch_dim.y;
61 float dx = 1.f / float(NX);
62 float dy = 1.f / float(NY);
63
64 // Map sample to rectangle [-0.5,0.5] [-0.5,0.5]
65 sp.x = -0.5f + ii * dx + float(launch_index.x) * dx / float(Nx) + Rx * dx / float(Nx);
66 sp.y = -0.5f + jj * dy + float(launch_index.y) * dy / float(Ny) + Ry * dy / float(Ny);
67 sp.z = 0.f;
68
69 int ID = maskID[objID];
70 if (ID >= 0 && primitive_solid_fraction[UUID] > 0.f && primitive_solid_fraction[UUID] < 1.f) { // has texture transparency
71
72 d_sampleTexture_patch(sp, optix::make_int2(ii, jj), optix::make_float2(dx, dy), prd, ID, puvID);
73 }
74
75 // calculate rectangle normal vector (world coordinates)
76 float3 v0 = make_float3(0, 0, 0);
77 d_transformPoint(m_trans, v0);
78 float3 v1 = make_float3(1, 0, 0);
79 d_transformPoint(m_trans, v1);
80 float3 v2 = make_float3(0, 1, 0);
81 d_transformPoint(m_trans, v2);
82
83 normal = normalize(cross(v1 - v0, v2 - v0));
84
85 } else if (ptype == 1) { // Triangle
86
87 // Map sample to triangle with vertices (0,0,0), (0,1,0), (1,1,0)
88 if (Rx < Ry) {
89 sp.x = Rx;
90 sp.y = Ry;
91 } else {
92 sp.x = Ry;
93 sp.y = Rx;
94 }
95 sp.z = 0;
96
97 // calculate triangle normal vector (world coordinates)
98 float3 v0 = make_float3(0, 0, 0);
99 d_transformPoint(m_trans, v0);
100 float3 v1 = make_float3(0, 1, 0);
101 d_transformPoint(m_trans, v1);
102 float3 v2 = make_float3(1, 1, 0);
103 d_transformPoint(m_trans, v2);
104
105 normal = normalize(cross(v1 - v0, v2 - v1));
106
107 int ID = maskID[objID];
108 if (ID >= 0 && primitive_solid_fraction[UUID] > 0.f && primitive_solid_fraction[UUID] < 1.f) { // has texture transparency
109
110 d_sampleTexture_triangle(sp, v0, v1, v2, prd, m_trans, ID, puvID);
111 }
112
113 } else if (ptype == 2) { // Disk
114
115 d_sampleDisk(prd.seed, sp);
116
117 // calculate disk normal vector (world coordinates)
118 float3 v0 = make_float3(0, 0, 0);
119 d_transformPoint(m_trans, v0);
120 float3 v1 = make_float3(1, 0, 0);
121 d_transformPoint(m_trans, v1);
122 float3 v2 = make_float3(0, 1, 0);
123 d_transformPoint(m_trans, v2);
124
125 normal = normalize(cross(v1 - v0, v2 - v0));
126
127 } else if (ptype == 4) { // Voxel
128
129 float Rz = rnd(prd.seed);
130 }
131
132 // translate the ray to the location of the primitive
133
134 float3 ray_origin = sp;
135 d_transformPoint(m_trans, ray_origin);
136
137 // Send a ray toward each source
138 for (int rr = 0; rr < Nsources; rr++) {
139
140 // set the ray direction
141 float3 ray_direction;
142 float ray_magnitude;
143 if (source_types[rr] == 0) { // collimated source
144 ray_direction = normalize(source_positions[rr]);
145 ray_magnitude = RT_DEFAULT_MAX;
146 prd.strength = 1. / double(launch_dim.x * launch_dim.y) * fabs(dot(normal, ray_direction));
147 } else if (source_types[rr] == 1 || source_types[rr] == 2) { // sphere source
148
149 // sample point on surface of sphere
150 float theta_s = acos_safe(1.f - 2.f * rnd(prd.seed));
151 float phi_s = rnd(prd.seed) * 2.f * M_PI;
152 float3 sphere_point = 0.5 * source_widths[rr].x * make_float3(sin(theta_s) * cos(phi_s), sin(theta_s) * sin(phi_s), cos(theta_s));
153
154 ray_direction = sphere_point + source_positions[rr] - ray_origin;
155
156 ray_magnitude = d_magnitude(ray_direction);
157 ray_direction = normalize(ray_direction);
158 prd.strength = 0.f;
159 uint N = 10;
160 for (uint j = 0; j < N; j++) {
161 for (uint i = 0; i < N; i++) {
162 float theta = acos_safe(1.f - 2.f * (float(i) + 0.5f) / float(N));
163 float phi = (float(j) + 0.5f) * 2.f * M_PI / float(N);
164 float3 light_direction = make_float3(sin(theta) * cos(phi), sin(theta) * sin(phi), cos(theta));
165 if (dot(light_direction, ray_direction) < 0) {
166 prd.strength +=
167 1. / double(launch_dim.x * launch_dim.y) * fabs(dot(normal, ray_direction)) * fabs(dot(light_direction, ray_direction)) / (ray_magnitude * ray_magnitude) / (N * N) * source_widths[rr].x * source_widths[rr].x;
168 }
169 }
170 }
171
172 } else if (source_types[rr] == 3) { // rectangle source
173
174 // transformation matrix
175 float light_transform[16];
176 d_makeTransformMatrix(source_rotations[rr], light_transform);
177
178 // sample point on surface of disk
179 float3 square_point;
180 d_sampleSquare(prd.seed, square_point);
181 square_point = make_float3(source_widths[rr].x * square_point.x, source_widths[rr].y * square_point.y, square_point.z);
182 d_transformPoint(light_transform, square_point);
183
184 float3 light_direction = make_float3(0, 0, 1);
185 d_transformPoint(light_transform, light_direction);
186
187 ray_direction = square_point + source_positions[rr] - ray_origin;
188
189 if (dot(ray_direction, light_direction) > 0.f) { // don't emit from back side of light source (note that ray goes toward the source, so the dot produce is negative when light is pointed at primitive)
190 continue;
191 }
192
193 ray_magnitude = d_magnitude(ray_direction);
194 ray_direction = normalize(ray_direction);
195 prd.strength = 1. / double(launch_dim.x * launch_dim.y) * fabs(dot(normal, ray_direction)) * fabs(dot(light_direction, ray_direction)) / (ray_magnitude * ray_magnitude) * source_widths[rr].x * source_widths[rr].y / M_PI;
196
197 } else if (source_types[rr] == 4) { // disk source
198
199 // transformation matrix
200 float light_transform[16];
201 d_makeTransformMatrix(source_rotations[rr], light_transform);
202
203 // sample point on surface of disk
204 float3 disk_point;
205 d_sampleDisk(prd.seed, disk_point);
206 d_transformPoint(light_transform, disk_point);
207
208 float3 light_direction = make_float3(0, 0, 1);
209 d_transformPoint(light_transform, light_direction);
210
211 ray_direction = source_widths[rr].x * disk_point + source_positions[rr] - ray_origin;
212
213 if (dot(ray_direction, light_direction) > 0.f) { // don't emit from back side of light source (note that ray goes toward the source, so the dot produce is negative when light is pointed at primitive)
214 continue;
215 }
216
217 ray_magnitude = d_magnitude(ray_direction);
218 ray_direction = normalize(ray_direction);
219 prd.strength = 1. / double(launch_dim.x * launch_dim.y) * fabs(dot(normal, ray_direction)) * fabs(dot(light_direction, ray_direction)) / (ray_magnitude * ray_magnitude) * source_widths[rr].x * source_widths[rr].x;
220 }
221
222 optix::Ray ray = optix::make_Ray(ray_origin, ray_direction, direct_ray_type, 1e-4, ray_magnitude);
223
224 prd.origin_UUID = UUID;
225 prd.source_ID = rr;
226 prd.hit_periodic_boundary = false;
227
228 if (dot(ray_direction, normal) > 0) {
229 prd.face = 1;
230 } else {
231 prd.face = 0;
232 }
233
234 if ((prd.face == 1 || twosided_flag[objID] == 1) && twosided_flag[objID] != 3) {
235
236 for (int wrap = 0; wrap < 10; ++wrap) {
237 rtTrace(top_object, ray, prd);
238
239 if (!prd.hit_periodic_boundary)
240 break; // real hit or miss → done
241
242 ray.origin = prd.periodic_hit;
243 prd.hit_periodic_boundary = false;
244 }
245 }
246 }
247 }
248 }
249}
250
251RT_PROGRAM void diffuse_raygen() {
252
253 uint dimx = launch_dim.x * launch_dim.y;
254 uint indx = launch_dim.x * launch_index.y + launch_index.x;
255
256 PerRayData prd;
257 prd.seed = tea<16>(indx + dimx * launch_index.z, random_seed);
258
259 uint objID = launch_offset + launch_index.z;
260
261 if (launch_face == 0 && twosided_flag[objID] == 0) { // skip the launch if from the bottom face and twosided_flag = 0
262 return;
263 }
264
265 uint pID = primitiveID[objID];
266 uint ptype = primitive_type[objID];
267 int puvID = uvID[objID];
268
269 // transformation matrix
270 float m_trans[16];
271 for (uint i = 0; i < 16; i++) {
272 m_trans[i] = transform_matrix[optix::make_uint2(i, objID)];
273 }
274
275 float3 sp, normal;
276
277 // looping over sub-patches
278 int NX = object_subdivisions[objID].x;
279 int NY = object_subdivisions[objID].y;
280 for (int jj = 0; jj < NY; jj++) {
281 for (int ii = 0; ii < NX; ii++) {
282
283 uint UUID = pID + jj * NX + ii;
284
285 // two random samples [0,1]
286 float Rx = rnd(prd.seed);
287 float Ry = rnd(prd.seed);
288
289 if (ptype > 4) {
290 printf("objID = %d\n", objID);
291 printf("Invalid primitive type in diffuse ray launch.\n");
292 }
293
294 if (ptype == 0 || ptype == 3) { // Patch or Tile
295
296 // calculate rectangle normal vector (world coordinates)
297 float3 s0 = make_float3(0, 0, 0);
298 float3 s1 = make_float3(1, 0, 0);
299 float3 s2 = make_float3(0, 1, 0);
300 d_transformPoint(m_trans, s0);
301 d_transformPoint(m_trans, s1);
302 d_transformPoint(m_trans, s2);
303
304 normal = normalize(cross(s1 - s0, s2 - s0));
305
306 float dx = 1.f / float(NX);
307 float dy = 1.f / float(NY);
308
309 // Map sample to rectangle [-0.5,0.5] [-0.5,0.5]
310 sp.x = -0.5f + (ii + Rx) * dx;
311 sp.y = -0.5f + (jj + Ry) * dy;
312 sp.z = 0.f;
313
314 int ID = maskID[objID];
315 if (ID >= 0) { // has texture transparency
316
317 d_sampleTexture_patch(sp, optix::make_int2(ii, jj), optix::make_float2(dx, dy), prd, ID, puvID);
318 }
319
320 } else if (ptype == 1) { // Triangle
321
322 // Map sample to triangle with vertices (0,0,0), (0,1,0), (1,1,0)
323 if (Rx < Ry) {
324 sp.x = Rx;
325 sp.y = Ry;
326 } else {
327 sp.x = Ry;
328 sp.y = Rx;
329 }
330 sp.z = 0;
331
332 // calculate triangle normal vector (world coordinates)
333 float3 v0 = make_float3(0, 0, 0);
334 d_transformPoint(m_trans, v0);
335 float3 v1 = make_float3(0, 1, 0);
336 d_transformPoint(m_trans, v1);
337 float3 v2 = make_float3(1, 1, 0);
338 d_transformPoint(m_trans, v2);
339
340 normal = normalize(cross(v1 - v0, v2 - v1));
341
342 int ID = maskID[objID];
343 if (ID >= 0) { // has texture transparency
344
345 d_sampleTexture_triangle(sp, v0, v1, v2, prd, m_trans, ID, puvID);
346 }
347
348 } else if (ptype == 2) { // Disk
349
350 // Map Sample to disk - from Suffern (2007) "Ray tracing fom the ground up" Chap. 6
351
352 // first map sample point to rectangle [-1,1] [-1,1]
353 sp.x = -1.f + 2.f * Rx;
354 sp.y = -1.f + 2.f * Ry;
355
356 float r, p;
357 if (sp.x > -sp.y) {
358 if (sp.x > sp.y) {
359 r = sp.x;
360 p = sp.y / sp.x;
361 } else {
362 r = sp.y;
363 p = 2.f - sp.x / sp.y;
364 }
365 } else {
366 if (sp.x < sp.y) {
367 r = -sp.x;
368 p = 4.f + sp.y / sp.x;
369 } else {
370 r = -sp.y;
371 if (sp.y != 0.f) { // avoid division by zero at origin
372 p = 6.f - sp.x / sp.y;
373 } else {
374 p = 0.f;
375 }
376 }
377 }
378 p *= 0.25f * M_PI;
379
380 // find x,y point on unit disk
381 sp.x = r * cosf(p);
382 sp.y = r * sinf(p);
383 sp.z = 0.f;
384
385 // calculate disk normal vector (world coordinates)
386 float3 v0 = make_float3(0, 0, 0);
387 d_transformPoint(m_trans, v0);
388 float3 v1 = make_float3(1, 0, 0);
389 d_transformPoint(m_trans, v1);
390 float3 v2 = make_float3(0, 1, 0);
391 d_transformPoint(m_trans, v2);
392 normal = normalize(cross(v1 - v0, v2 - v0));
393
394 } else if (ptype == 4) { // Voxel
395
396 // Map sample to cube [-0.5,0.5] [-0.5,0.5] [-0.5,0.5]
397 sp.x = -0.5f + Rx;
398 sp.y = -0.5f + Ry;
399 sp.z = -0.5f + rnd(prd.seed);
400 }
401
402 // Choose random hemispherical direction - map samples to hemisphere (from Suffern (2007) "Ray tracing fom the ground up" Chap. 6)
403
404 float Rt;
405 float Rp;
406
407 Rt = (launch_index.x + rnd(prd.seed)) / float(launch_dim.x);
408 Rp = (launch_index.y + rnd(prd.seed)) / float(launch_dim.y);
409
410 float t;
411 if (ptype == 4) { // voxel
412 t = acos_safe(1.f - Rt);
413 } else { // other
414 t = asin_safe(sqrtf(Rt));
415 }
416 float p = 2.f * M_PI * Rp;
417
418 float3 ray_direction;
419 ray_direction.x = sin(t) * cos(p);
420 ray_direction.y = sin(t) * sin(p);
421 ray_direction.z = cos(t);
422
423 float3 ray_origin;
424 optix::Ray ray;
425
426 if (ptype == 4) { // voxel
427
428 prd.strength = 0.5f / float(dimx);
429 prd.origin_UUID = UUID;
430 prd.face = 0;
431 prd.source_ID = 0;
432 prd.hit_periodic_boundary = false;
433
434 ray_origin = sp;
435 d_transformPoint(m_trans, ray_origin);
436
437 ray = optix::make_Ray(ray_origin, ray_direction, diffuse_ray_type, 1e-5, RT_DEFAULT_MAX);
438 rtTrace(top_object, ray, prd);
439
440 ray = optix::make_Ray(ray_origin, -ray_direction, diffuse_ray_type, 1e-5, RT_DEFAULT_MAX);
441 rtTrace(top_object, ray, prd);
442
443 } else { // not a voxel
444
445 ray_direction = d_rotatePoint(ray_direction, acos_safe(normal.z), atan2(normal.y, normal.x));
446
447 prd.strength = 1.f / float(dimx);
448
449 prd.origin_UUID = UUID;
450 prd.source_ID = 0;
451 prd.hit_periodic_boundary = false;
452
453 // ---- "top" surface launch -------
454 ray_origin = sp;
455 d_transformPoint(m_trans, ray_origin);
456
457 if (launch_face == 1 && twosided_flag[objID] != 3) {
458
459 ray = optix::make_Ray(ray_origin, ray_direction, diffuse_ray_type, 1e-5, RT_DEFAULT_MAX);
460
461 prd.face = 1;
462
463 for (int wrap = 0; wrap < 10; ++wrap) {
464 rtTrace(top_object, ray, prd);
465
466 if (!prd.hit_periodic_boundary)
467 break; // real hit or miss → done
468
469 ray.origin = prd.periodic_hit;
470 prd.hit_periodic_boundary = false;
471 }
472
473 // ---- "bottom" surface launch -------
474 } else if (launch_face == 0 && twosided_flag[objID] == 1) {
475 ray_direction = -ray_direction;
476 ray = optix::make_Ray(ray_origin, ray_direction, diffuse_ray_type, 1e-5, RT_DEFAULT_MAX);
477
478 prd.face = 0;
479
480 for (int wrap = 0; wrap < 10; ++wrap) {
481 rtTrace(top_object, ray, prd);
482
483 if (!prd.hit_periodic_boundary)
484 break; // real hit or miss → done
485
486 ray.origin = prd.periodic_hit;
487 prd.hit_periodic_boundary = false;
488 }
489 }
490 }
491 }
492 }
493}
494
495RT_PROGRAM void camera_raygen() {
496
497 uint dimx = launch_dim.x * launch_dim.y; // x number of ray, y width, z length
498 uint indx = launch_dim.x * launch_index.y + launch_index.x;
499
500 optix::int2 camera_resolution = optix::make_int2(launch_dim.y, launch_dim.z);
501
502 PerRayData prd;
503 prd.seed = tea<16>(indx + dimx * launch_index.z, random_seed);
504
505 float3 sp;
506
507 uint ii = launch_index.y; // x-pixel
508 uint jj = launch_index.z; // y-pixel
509 size_t origin_ID = jj * launch_dim.y + ii; // global pixel index
510
511
512 // distortion
513 // float PPointsRatiox =1.052f;
514 // float PPointsRatioy =0.999f;
515 // float sensorxscale = 1.0054;
516 // float focalxy = 710;
517 // double x =(float(ii)-camera_resolution.x/2 * PPointsRatiox)/focalxy*sensorxscale;// / focalxy; cam_res.y = 712
518 // double y = (float(jj)-camera_resolution.y/2 * PPointsRatioy )/focalxy; /// focalxy; cam_res.x = 1072
519 // double r2 = x*x + y*y;
520 // double distCoeffs[4] = {-0.3535674,0.17298, 0, 0};
521 // double ii_d = x * (1+ distCoeffs[0] * r2 + distCoeffs[1] * r2 * r2) + 2 * distCoeffs[2] * x * y + distCoeffs[3] * (r2 + 2 * x * x);
522 // double jj_d = y * (1+ distCoeffs[0] * r2 + distCoeffs[1] * r2 * r2) + 2 * distCoeffs[3] * x * y + distCoeffs[2] * (r2 + 2 * y * y);
523 // ii_d = ii_d*focalxy+float(camera_resolution.x)/2 * PPointsRatiox;
524 // jj_d = jj_d*focalxy+float(camera_resolution.y)/2 * PPointsRatioy;
525
526 // *** sample a point on the pixel (view direction coordinate aligned) *** //
527
528 float PPointsRatiox = 1.f;
529 float PPointsRatioy = 1.f;
530 float Rx = rnd(prd.seed);
531 float Ry = rnd(prd.seed);
532
533 // Map sample to pixel
534 // Calculate VFOV scaling factor directly from aspect ratio
535 // Since FOV_aspect_ratio = HFOV/VFOV, then tan(half_VFOV) = tan(half_HFOV)/FOV_aspect_ratio
536 // The multiplier = tan(half_VFOV)/tan(half_HFOV) = 1/FOV_aspect_ratio
537 float multiplier = 1.0f / FOV_aspect_ratio;
538 sp.y = (-0.5f * PPointsRatioy + (ii + Rx) / float(camera_resolution.x));
539 sp.z = (0.5f * PPointsRatiox - (jj + Ry) / float(camera_resolution.y)) * multiplier;
540 sp.x = camera_viewplane_length;
541
542 // *** Determine point 'p' on focal plane that passes through the lens center (0,0) and pixel sample (view direction coordinate aligned) *** //
543
544 float3 p = make_float3(camera_focal_length, sp.y / camera_viewplane_length * camera_focal_length, sp.z / camera_viewplane_length * camera_focal_length);
545
546 // *** Sample point on lens (view direction coordinate aligned) *** //
547
548 float3 ray_origin = make_float3(0, 0, 0);
549 if (camera_lens_diameter > 0) {
550 float3 disk_sample;
551 d_sampleDisk(prd.seed, disk_sample);
552 ray_origin = make_float3(0.f, 0.5f * disk_sample.x * camera_lens_diameter, 0.5f * disk_sample.y * camera_lens_diameter);
553 }
554
555 //*** ray direction is line from lens sample to p ***//
556
557 float3 ray_direction = p - ray_origin;
558
559 //*** rotate ray origin and direction into the direction of the camera view *** //
560
561 ray_origin = d_rotatePoint(ray_origin, -0.5 * M_PI + camera_direction.x, 0.5f * M_PI - camera_direction.y) + camera_position;
562
563 ray_direction = d_rotatePoint(ray_direction, -0.5 * M_PI + camera_direction.x, 0.5f * M_PI - camera_direction.y);
564 ray_direction /= d_magnitude(ray_direction);
565
566 optix::Ray ray;
567
568 prd.strength = 1.f / float(launch_dim.x);
569
570 prd.origin_UUID = origin_ID;
571 prd.face = 1;
572 prd.source_ID = 0;
573 prd.hit_periodic_boundary = false;
574
575 ray = optix::make_Ray(ray_origin, ray_direction, camera_ray_type, 1e-5, RT_DEFAULT_MAX);
576
577 for (int wrap = 0; wrap < 10; ++wrap) {
578 rtTrace(top_object, ray, prd);
579
580 if (!prd.hit_periodic_boundary)
581 break; // real hit or miss → done
582
583 ray.origin = prd.periodic_hit;
584 prd.hit_periodic_boundary = false;
585 }
586}
587
588RT_PROGRAM void pixel_label_raygen() {
589
590 uint indx = launch_dim.y * launch_index.z + launch_index.y;
591
592 optix::int2 camera_resolution = optix::make_int2(launch_dim.y, launch_dim.z);
593
594 PerRayData prd;
595 prd.seed = tea<16>(indx, random_seed);
596
597 float3 sp;
598
599 uint ii = launch_index.y; // x-pixel
600
601 uint jj = launch_index.z; // y-pixel
602
603 size_t origin_ID = jj * launch_dim.y + ii; // global pixel index
604
605 // Map sample to center of pixel
606 // Calculate VFOV scaling factor directly from aspect ratio
607 // Since FOV_aspect_ratio = HFOV/VFOV, then tan(half_VFOV) = tan(half_HFOV)/FOV_aspect_ratio
608 // The multiplier = tan(half_VFOV)/tan(half_HFOV) = 1/FOV_aspect_ratio
609 float multiplier = 1.0f / FOV_aspect_ratio;
610 sp.y = (-0.5f + (ii + 0.5f) / float(camera_resolution.x));
611 sp.z = (0.5f - (jj + 0.5f) / float(camera_resolution.y)) * multiplier;
612 sp.x = camera_viewplane_length;
613
614
615 // *** Determine point 'p' on focal plane that passes through the lens center (0,0) and pixel sample (view direction coordinate aligned) *** //
616
617 float3 p = make_float3(camera_focal_length, sp.y / camera_viewplane_length * camera_focal_length, sp.z / camera_viewplane_length * camera_focal_length);
618
619 // *** Ray is launched from center of lens *** //
620
621 float3 ray_origin = make_float3(0, 0, 0);
622
623 //*** ray direction is line from ray origin to p ***//
624
625 float3 ray_direction = p;
626
627 //*** rotate ray origin and direction into the direction of the camera view *** //
628
629 ray_origin = d_rotatePoint(ray_origin, -0.5 * M_PI + camera_direction.x, 0.5f * M_PI - camera_direction.y) + camera_position;
630
631 ray_direction = d_rotatePoint(ray_direction, -0.5 * M_PI + camera_direction.x, 0.5f * M_PI - camera_direction.y);
632 ray_direction /= d_magnitude(ray_direction);
633
634 optix::Ray ray;
635
636 prd.strength = 0.f;
637
638 prd.origin_UUID = origin_ID;
639 prd.face = 1;
640 prd.source_ID = 0;
641 prd.hit_periodic_boundary = false;
642
643 ray = optix::make_Ray(ray_origin, ray_direction, pixel_label_ray_type, 1e-5, RT_DEFAULT_MAX);
644
645 for (int wrap = 0; wrap < 10; ++wrap) {
646 rtTrace(top_object, ray, prd);
647
648 if (!prd.hit_periodic_boundary)
649 break; // real hit or miss → done
650
651 ray.origin = prd.periodic_hit;
652 prd.hit_periodic_boundary = false;
653 }
654}