sdis_solve_camera.c (22497B)
1 /* Copyright (C) 2016-2025 |Méso|Star> (contact@meso-star.com) 2 * 3 * This program is free software: you can redistribute it and/or modify 4 * it under the terms of the GNU General Public License as published by 5 * the Free Software Foundation, either version 3 of the License, or 6 * (at your option) any later version. 7 * 8 * This program is distributed in the hope that it will be useful, 9 * but WITHOUT ANY WARRANTY; without even the implied warranty of 10 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 11 * GNU General Public License for more details. 12 * 13 * You should have received a copy of the GNU General Public License 14 * along with this program. If not, see <http://www.gnu.org/licenses/>. */ 15 16 #include "sdis.h" 17 #include "sdis_c.h" 18 #include "sdis_camera.h" 19 #include "sdis_device_c.h" 20 #include "sdis_estimator_buffer_c.h" 21 #include "sdis_log.h" 22 #include "sdis_medium_c.h" 23 #include "sdis_realisation.h" 24 #include "sdis_scene_c.h" 25 #include "sdis_tile.h" 26 #ifdef SDIS_ENABLE_MPI 27 #include "sdis_mpi.h" 28 #endif 29 30 #include <rsys/clock_time.h> 31 #include <rsys/cstr.h> 32 #include <rsys/list.h> 33 #include <rsys/morton.h> 34 35 #include <omp.h> 36 37 /******************************************************************************* 38 * Helper function 39 ******************************************************************************/ 40 static res_T 41 check_solve_camera_args(const struct sdis_solve_camera_args* args) 42 { 43 if(!args) return RES_BAD_ARG; 44 if(!args->cam) return RES_BAD_ARG; 45 46 /* Check the image resolution */ 47 if(!args->image_definition[0] || !args->image_definition[1]) { 48 return RES_BAD_ARG; 49 } 50 51 /* Check the number of samples per pixel */ 52 if(!args->spp) { 53 return RES_BAD_ARG; 54 } 55 56 /* Check the time range */ 57 if(args->time_range[0] < 0 || args->time_range[1] < args->time_range[0]) { 58 return RES_BAD_ARG; 59 } 60 if(args->time_range[1] > DBL_MAX 61 && args->time_range[0] != args->time_range[1]) { 62 return RES_BAD_ARG; 63 } 64 65 /* Check the picard order */ 66 if(args->picard_order < 1) { 67 return RES_BAD_ARG; 68 } 69 70 /* Check RNG type */ 71 if(!args->rng_state && args->rng_type >= SSP_RNG_TYPES_COUNT__) { 72 return RES_BAD_ARG; 73 } 74 75 /* Check the diffusion algorithm */ 76 if((unsigned)args->diff_algo >= SDIS_DIFFUSION_ALGORITHMS_COUNT__) { 77 return RES_BAD_ARG; 78 } 79 80 return RES_OK; 81 } 82 83 static res_T 84 solve_pixel 85 (struct sdis_scene* scn, 86 struct ssp_rng* rng, 87 const unsigned enc_id, 88 const struct sdis_camera* cam, 89 const double time_range[2], /* Observation time */ 90 const size_t ipix[2], /* Pixel coordinate in the image space */ 91 const size_t nrealisations, 92 const int register_paths, /* Combination of enum sdis_heat_path_flag */ 93 const double pix_sz[2], /* Pixel size in the normalized image plane */ 94 const size_t picard_order, 95 const enum sdis_diffusion_algorithm diff_algo, 96 struct sdis_estimator* estimator, 97 struct pixel* pixel) 98 { 99 struct sdis_heat_path* pheat_path = NULL; 100 size_t irealisation; 101 res_T res = RES_OK; 102 ASSERT(scn && rng && cam && ipix && nrealisations); 103 ASSERT(pix_sz && pix_sz[0] > 0 && pix_sz[1] > 0); 104 ASSERT(pixel && time_range); 105 106 pixel->acc_temp = ACCUM_NULL; 107 pixel->acc_time = ACCUM_NULL; 108 109 FOR_EACH(irealisation, 0, nrealisations) { 110 struct ray_realisation_args realis_args = RAY_REALISATION_ARGS_NULL; 111 struct time t0, t1; 112 double samp[2]; /* Pixel sample */ 113 double ray_pos[3]; 114 double ray_dir[3]; 115 double w = 0; 116 struct sdis_heat_path heat_path; 117 double time; 118 res_T res_simul = RES_OK; 119 120 /* Begin time registration */ 121 time_current(&t0); 122 123 time = sample_time(rng, time_range); 124 if(register_paths) { 125 heat_path_init(scn->dev->allocator, &heat_path); 126 pheat_path = &heat_path; 127 } 128 129 /* Generate a sample into the pixel to estimate */ 130 samp[0] = ((double)ipix[0] + ssp_rng_canonical(rng)) * pix_sz[0]; 131 samp[1] = ((double)ipix[1] + ssp_rng_canonical(rng)) * pix_sz[1]; 132 133 /* Generate a ray starting from the camera position and passing through 134 * pixel sample */ 135 camera_ray(cam, samp, ray_pos, ray_dir); 136 137 /* Launch the realisation */ 138 realis_args.rng = rng; 139 realis_args.enc_id = enc_id; 140 realis_args.time = time; 141 realis_args.picard_order = picard_order; 142 realis_args.heat_path = pheat_path; 143 realis_args.irealisation = (size_t)irealisation; 144 realis_args.diff_algo = diff_algo; 145 d3_set(realis_args.position, ray_pos); 146 d3_set(realis_args.direction, ray_dir); 147 res_simul = ray_realisation_3d(scn, &realis_args, &w); 148 149 /* Handle fatal error */ 150 if(res_simul != RES_OK && res_simul != RES_BAD_OP) { 151 res = res_simul; 152 goto error; 153 } 154 155 if(pheat_path) { 156 pheat_path->status = res_simul == RES_OK 157 ? SDIS_HEAT_PATH_SUCCESS 158 : SDIS_HEAT_PATH_FAILURE; 159 160 /* Check if the path must be saved regarding the register_paths mask */ 161 if(!(register_paths & (int)pheat_path->status)) { 162 heat_path_release(pheat_path); 163 pheat_path = NULL; 164 } else { /* Register the sampled path */ 165 ASSERT(estimator); 166 res = estimator_add_and_release_heat_path(estimator, pheat_path); 167 if(res != RES_OK) goto error; 168 pheat_path = NULL; 169 } 170 } 171 172 /* Stop time registration */ 173 time_sub(&t0, time_current(&t1), &t0); 174 175 if(res_simul == RES_OK) { 176 /* Update pixel accumulators */ 177 const double usec = (double)time_val(&t0, TIME_NSEC) * 0.001; 178 pixel->acc_temp.sum += w; 179 pixel->acc_temp.sum2 += w*w; 180 pixel->acc_temp.count += 1; 181 pixel->acc_time.sum += usec; 182 pixel->acc_time.sum2 += usec*usec; 183 pixel->acc_time.count += 1; 184 } 185 } 186 187 exit: 188 if(pheat_path) heat_path_release(pheat_path); 189 return res; 190 error: 191 goto exit; 192 } 193 194 static res_T 195 solve_tile 196 (struct sdis_scene* scn, 197 struct ssp_rng* rng, 198 const unsigned enc_id, 199 const struct sdis_camera* cam, 200 const double time_range[2], 201 const size_t tile_org[2], /* Origin of the tile in pixel space */ 202 const size_t tile_size[2], /* #pixels in the tile in X and Y */ 203 const size_t spp, /* #samples per pixel */ 204 const int register_paths, /* Combination of enum sdis_heat_path_flag */ 205 const double pix_sz[2], /* Pixel size in the normalized image plane */ 206 const size_t picard_order, 207 const enum sdis_diffusion_algorithm diff_algo, 208 struct sdis_estimator_buffer* buf, 209 struct tile* tile) 210 { 211 size_t mcode; /* Morton code of the tile pixel */ 212 size_t npixels; 213 res_T res = RES_OK; 214 ASSERT(scn && rng && cam && spp); 215 ASSERT(tile_size && tile_size[0] && tile_size[1]); 216 ASSERT(pix_sz && pix_sz[0] > 0 && pix_sz[1] > 0 && time_range); 217 218 /* Adjust the #pixels to process them wrt a morton order */ 219 npixels = round_up_pow2(MMAX(tile_size[0], tile_size[1])); 220 npixels *= npixels; 221 222 FOR_EACH(mcode, 0, npixels) { 223 struct pixel* pixel = NULL; 224 struct sdis_estimator* estimator = NULL; 225 uint16_t ipix_tile[2]; 226 size_t ipix_image[2]; 227 228 ipix_tile[0] = morton2D_decode_u16((uint32_t)(mcode>>0)); 229 if(ipix_tile[0] >= tile_size[0]) continue; 230 ipix_tile[1] = morton2D_decode_u16((uint32_t)(mcode>>1)); 231 if(ipix_tile[1] >= tile_size[1]) continue; 232 233 pixel = tile_at(tile, ipix_tile[0], ipix_tile[1]); 234 235 /* Compute the pixel coordinates in image space */ 236 ipix_image[0] = ipix_tile[0] + tile_org[0]; 237 ipix_image[1] = ipix_tile[1] + tile_org[1]; 238 239 if(register_paths != SDIS_HEAT_PATH_NONE) { 240 ASSERT(buf); 241 estimator = estimator_buffer_grab(buf, ipix_image[0], ipix_image[1]); 242 } 243 res = solve_pixel 244 (scn, rng, enc_id, cam, time_range, ipix_image, spp, register_paths, 245 pix_sz, picard_order, diff_algo, estimator, pixel); 246 if(res != RES_OK) goto error; 247 } 248 249 exit: 250 return res; 251 error: 252 goto exit; 253 } 254 255 static INLINE res_T 256 setup_estimator_from_pixel 257 (struct sdis_estimator* estimator, 258 const size_t spp, /* #samples per pixel */ 259 struct pixel* pixel) 260 { 261 ASSERT(estimator && spp && pixel); 262 ASSERT(pixel->acc_temp.count == pixel->acc_time.count); 263 estimator_setup_realisations_count 264 (estimator, spp, pixel->acc_temp.count); 265 estimator_setup_temperature 266 (estimator, pixel->acc_temp.sum, pixel->acc_temp.sum2); 267 estimator_setup_realisation_time 268 (estimator, pixel->acc_time.sum, pixel->acc_time.sum2); 269 return RES_OK; 270 } 271 272 static res_T 273 write_tile 274 (struct sdis_estimator_buffer* buf, 275 const size_t spp, /* #samples per pixel */ 276 struct tile* tile) 277 { 278 res_T res = RES_OK; 279 size_t tile_org[2]; 280 size_t buf_sz[2]; 281 size_t tile_sz[2]; 282 uint16_t x, y; 283 ASSERT(buf && spp && tile); 284 285 SDIS(estimator_buffer_get_definition(buf, buf_sz)); 286 287 tile_org[0] = (size_t)(tile->data.x * TILE_SIZE); 288 tile_org[1] = (size_t)(tile->data.y * TILE_SIZE); 289 tile_sz[0] = MMIN(TILE_SIZE, buf_sz[0] - tile_org[0]); 290 tile_sz[1] = MMIN(TILE_SIZE, buf_sz[1] - tile_org[1]); 291 292 FOR_EACH(y, 0, tile_sz[1]) { 293 const size_t pix_y = tile_org[1] + y; 294 FOR_EACH(x, 0, tile_sz[0]) { 295 const size_t pix_x = tile_org[0] + x; 296 struct sdis_estimator* estimator = NULL; 297 struct pixel* pixel = NULL; 298 299 estimator = estimator_buffer_grab(buf, pix_x, pix_y); 300 pixel = tile_at(tile, x, y); 301 302 res = setup_estimator_from_pixel(estimator, spp, pixel); 303 if(res != RES_OK) goto error; 304 } 305 } 306 307 exit: 308 return res; 309 error: 310 goto exit; 311 } 312 313 static res_T 314 write_list_of_tiles 315 (struct sdis_estimator_buffer* buf, 316 const size_t spp, /* #samples per pixel */ 317 struct list_node* tiles) /* Tiles to write */ 318 { 319 struct list_node* node = NULL; 320 res_T res = RES_OK; 321 ASSERT(buf && spp && tiles); 322 323 LIST_FOR_EACH(node, tiles) { 324 struct tile* tile = CONTAINER_OF(node, struct tile, node); 325 res = write_tile(buf, spp, tile); 326 if(res != RES_OK) goto error; 327 } 328 329 exit: 330 return res; 331 error: 332 goto exit; 333 } 334 335 #ifndef SDIS_ENABLE_MPI 336 static INLINE res_T 337 gather_tiles 338 (struct sdis_device* dev, 339 struct sdis_estimator_buffer* buf, /* NULL on non master processes */ 340 const size_t spp, /* #realisations per pixel */ 341 const size_t ntiles, /* Overall #tiles that must be written into `buf' */ 342 struct list_node* tiles) /* List of tiles of the current process */ 343 { 344 (void)dev, (void)ntiles; 345 return write_list_of_tiles(buf, spp, tiles); 346 } 347 #else 348 /* Gather the tiles and write them into the estimator buffer. The master process 349 * write its own tiles into the estimator buffer. Then, it gathers the tiles 350 * send by non master processes and write them too into the estimator buffer */ 351 static res_T 352 gather_tiles 353 (struct sdis_device* dev, 354 struct sdis_estimator_buffer* buf, /* NULL on non master processes */ 355 const size_t spp, /* #realisations per pixel */ 356 const size_t ntiles, /* Overall #tiles that must be written into `buf' */ 357 struct list_node* tiles) /* List of tiles of the current process */ 358 { 359 struct tile* tile_temp = NULL; 360 struct list_node* node = NULL; 361 res_T res = RES_OK; 362 ASSERT(dev && spp && tiles); 363 364 if(!dev->use_mpi) { 365 res = write_list_of_tiles(buf, spp, tiles); 366 if(res != RES_OK) goto error; 367 goto exit; /* No more to do */ 368 } 369 370 /* Non master process */ 371 if(dev->mpi_rank != 0) { 372 373 /* Send to the master process the list of tiles solved by the current 374 * process */ 375 LIST_FOR_EACH(node, tiles) { 376 struct tile* tile = CONTAINER_OF(node, struct tile, node); 377 mutex_lock(dev->mpi_mutex); 378 MPI(Send(&tile->data, sizeof(tile->data), MPI_CHAR, 0, MPI_SDIS_MSG_TILE, 379 MPI_COMM_WORLD)); 380 mutex_unlock(dev->mpi_mutex); 381 } 382 383 /* Master process */ 384 } else { 385 size_t itile; 386 size_t ntiles_master; 387 388 /* Write into the buffer the tiles solved by the master process itself */ 389 res = write_list_of_tiles(buf, spp, tiles); 390 if(res != RES_OK) goto error; 391 392 /* Create a temporary tile to store the tile sent by the non master 393 * processes */ 394 res = tile_create(dev->allocator, &tile_temp); 395 if(res != RES_OK) { 396 log_err(dev, 397 "Could not allocate the tile to temporary store the tiles sent by " 398 "the non master processes -- %s", res_to_cstr(res)); 399 goto error; 400 } 401 402 /* Count the number of tiles rendered onto the master process */ 403 ntiles_master = 0; 404 LIST_FOR_EACH(node, tiles) ++ntiles_master; 405 ASSERT(ntiles_master <= ntiles); 406 407 /* Receive the remaining tiles sent by the non master processes */ 408 FOR_EACH(itile, ntiles_master, ntiles) { 409 MPI_Request req; 410 411 /* Asynchronously receive a tile */ 412 mutex_lock(dev->mpi_mutex); 413 MPI(Irecv(&tile_temp->data, sizeof(tile_temp->data), MPI_CHAR, 414 MPI_ANY_SOURCE, MPI_SDIS_MSG_TILE, MPI_COMM_WORLD, &req)); 415 mutex_unlock(dev->mpi_mutex); 416 mpi_waiting_for_request(dev, &req); 417 418 /* Write the tile into the estimator buffer */ 419 res = write_tile(buf, spp, tile_temp); 420 if(res != RES_OK) goto error; 421 } 422 } 423 424 exit: 425 if(tile_temp) tile_ref_put(tile_temp); 426 return res; 427 error: 428 goto exit; 429 } 430 #endif 431 432 /* Setup the accumulators of the whole estimator buffer */ 433 static res_T 434 finalize_estimator_buffer 435 (struct sdis_estimator_buffer* buf, 436 struct ssp_rng_proxy* rng_proxy, 437 const size_t spp) /* #samples per pixel */ 438 { 439 struct accum acc_temp = ACCUM_NULL; 440 struct accum acc_time = ACCUM_NULL; 441 size_t definition[2]; 442 size_t x, y; 443 size_t nrealisations = 0; 444 size_t nsuccesses = 0; 445 res_T res = RES_OK; 446 ASSERT(buf && rng_proxy && spp); 447 448 SDIS(estimator_buffer_get_definition(buf, definition)); 449 450 FOR_EACH(y, 0, definition[1]) { 451 FOR_EACH(x, 0, definition[0]) { 452 const struct sdis_estimator* estimator; 453 SDIS(estimator_buffer_at(buf, x, y, &estimator)); 454 acc_temp.sum += estimator->temperature.sum; 455 acc_temp.sum2 += estimator->temperature.sum2; 456 acc_temp.count += estimator->temperature.count; 457 acc_time.sum += estimator->realisation_time.sum; 458 acc_time.sum2 += estimator->realisation_time.sum2; 459 acc_time.count += estimator->realisation_time.count; 460 nsuccesses += estimator->nrealisations; 461 } 462 } 463 464 nrealisations = definition[0]*definition[1]*spp; 465 ASSERT(acc_temp.count == acc_time.count); 466 ASSERT(acc_temp.count == nsuccesses); 467 estimator_buffer_setup_realisations_count(buf, nrealisations, nsuccesses); 468 estimator_buffer_setup_temperature(buf, acc_temp.sum, acc_temp.sum2); 469 estimator_buffer_setup_realisation_time(buf, acc_time.sum, acc_time.sum2); 470 res = estimator_buffer_save_rng_state(buf, rng_proxy); 471 if(res != RES_OK) goto error; 472 473 exit: 474 return res; 475 error: 476 goto exit; 477 } 478 479 static void 480 free_rendered_tiles(struct list_node* tiles) 481 { 482 struct list_node* node; 483 struct list_node* tmp; 484 ASSERT(tiles); 485 LIST_FOR_EACH_SAFE(node, tmp, tiles) { 486 struct tile* tile = CONTAINER_OF(node, struct tile, node); 487 list_del(node); 488 tile_ref_put(tile); 489 } 490 } 491 492 /******************************************************************************* 493 * Exported function 494 ******************************************************************************/ 495 res_T 496 sdis_solve_camera 497 (struct sdis_scene* scn, 498 const struct sdis_solve_camera_args* args, 499 struct sdis_estimator_buffer** out_buf) 500 { 501 /* Time registration */ 502 struct time time0, time1; 503 char buffer[128]; /* Temporary buffer used to store formated time */ 504 505 /* Stardis variables */ 506 struct sdis_estimator_buffer* buf = NULL; 507 508 /* Random number generators */ 509 struct ssp_rng_proxy* rng_proxy = NULL; 510 struct ssp_rng** per_thread_rng = NULL; 511 512 /* Enclosure & medium in which the probe lies */ 513 unsigned enc_id = ENCLOSURE_ID_NULL; 514 515 /* Miscellaneous */ 516 size_t ntiles_x, ntiles_y, ntiles, ntiles_adjusted; 517 size_t ntiles_proc; /* #tiles for the current proc */ 518 struct list_node tiles; /* List of tiles rendered by the process */ 519 double pix_sz[2]; /* Size of a pixel in the normalized image plane */ 520 int64_t mcode; /* Morton code of a tile */ 521 int64_t mcode_1st; /* morton code of the 1st tile computed by the process */ 522 int64_t mcode_incr; /* Increment toward the next morton code */ 523 int32_t* progress = NULL; /* Per process progress bar */ 524 int pcent_progress = 1; /* Percentage requiring progress update */ 525 int register_paths = SDIS_HEAT_PATH_NONE; 526 int is_master_process = 1; 527 ATOMIC nsolved_tiles = 0; 528 ATOMIC res = RES_OK; 529 530 list_init(&tiles); 531 532 if(!scn || !out_buf) { res = RES_BAD_ARG; goto error; } 533 res = check_solve_camera_args(args); 534 if(res != RES_OK) goto error; 535 536 if(scene_is_2d(scn)) { 537 log_err(scn->dev, "%s: 2D scenes are not supported.\n", FUNC_NAME); 538 goto error; 539 } 540 541 /* Retrieve the medium in which the submitted position lies */ 542 res = scene_get_enclosure_id(scn, args->cam->position, &enc_id); 543 if(res != RES_OK) goto error; 544 545 /* Create the per thread RNGs */ 546 res = create_per_thread_rng 547 (scn->dev, args->rng_state, args->rng_type, &rng_proxy, &per_thread_rng); 548 if(res != RES_OK) goto error; 549 550 /* Allocate the per process progress status */ 551 res = alloc_process_progress(scn->dev, &progress); 552 if(res != RES_OK) goto error; 553 554 /* Compute the overall number of tiles */ 555 ntiles_x = (args->image_definition[0] + (TILE_SIZE-1)/*ceil*/)/TILE_SIZE; 556 ntiles_y = (args->image_definition[1] + (TILE_SIZE-1)/*ceil*/)/TILE_SIZE; 557 ntiles = ntiles_x * ntiles_y; 558 ntiles_adjusted = round_up_pow2(MMAX(ntiles_x, ntiles_y)); 559 ntiles_adjusted *= ntiles_adjusted; 560 561 #ifdef SDIS_ENABLE_MPI 562 is_master_process = !scn->dev->use_mpi || scn->dev->mpi_rank == 0; 563 if(scn->dev->use_mpi) { 564 mcode_1st = scn->dev->mpi_rank; 565 mcode_incr = scn->dev->mpi_nprocs; 566 567 /* Compute the #tiles of the current proc */ 568 ntiles_proc = ntiles / (size_t)scn->dev->mpi_nprocs; 569 if(ntiles % (size_t)scn->dev->mpi_nprocs > (size_t)scn->dev->mpi_rank) { 570 ++ntiles_proc; 571 } 572 } else 573 #endif 574 { 575 is_master_process = 1; 576 mcode_1st = 0; 577 mcode_incr = 1; 578 ntiles_proc = ntiles; 579 } 580 581 /* Update the progress bar every percent if escape sequences are allowed in 582 * log messages or only every 10 percent when only plain text is allowed. 583 * This reduces the number of lines of plain text printed */ 584 pcent_progress = scn->dev->no_escape_sequence ? 10 : 1; 585 586 /* Compute the normalized pixel size */ 587 pix_sz[0] = 1.0 / (double)args->image_definition[0]; 588 pix_sz[1] = 1.0 / (double)args->image_definition[1]; 589 590 /* Create the global estimator on the master process only */ 591 if(is_master_process) { 592 res = estimator_buffer_create 593 (scn->dev, args->image_definition[0], args->image_definition[1], &buf); 594 if(res != RES_OK) goto error; 595 } 596 597 /* Synchronise the processes */ 598 process_barrier(scn->dev); 599 600 #define PROGRESS_MSG "Rendering: " 601 print_progress(scn->dev, progress, PROGRESS_MSG); 602 603 /* Begin time registration of the computation */ 604 time_current(&time0); 605 606 /* Here we go! Launch the Monte Carlo estimation */ 607 omp_set_num_threads((int)scn->dev->nthreads); 608 register_paths = is_master_process 609 ? args->register_paths : SDIS_HEAT_PATH_NONE; 610 #pragma omp parallel for schedule(static, 1/*chunk size*/) 611 for(mcode = mcode_1st; mcode < (int64_t)ntiles_adjusted; mcode+=mcode_incr) { 612 size_t tile_org[2] = {0, 0}; 613 size_t tile_sz[2] = {0, 0}; 614 struct tile* tile = NULL; 615 const int ithread = omp_get_thread_num(); 616 struct ssp_rng* rng = per_thread_rng[ithread]; 617 size_t n; 618 int pcent; 619 res_T res_local = RES_OK; 620 621 if(ATOMIC_GET(&res) != RES_OK) continue; 622 623 tile_org[0] = morton2D_decode_u16((uint32_t)(mcode>>0)); 624 if(tile_org[0] >= ntiles_x) continue; /* Discard tile */ 625 tile_org[1] = morton2D_decode_u16((uint32_t)(mcode>>1)); 626 if(tile_org[1] >= ntiles_y) continue; /* Discard tile */ 627 628 res_local = tile_create(scn->dev->allocator, &tile); 629 if(tile == NULL) { 630 log_err(scn->dev, "%s: error allocating the tile (%lu, %lu) -- %s\n", 631 FUNC_NAME, 632 (unsigned long)tile_org[0], 633 (unsigned long)tile_org[1], 634 res_to_cstr(res_local)); 635 ATOMIC_SET(&res, res_local); 636 continue; 637 } 638 639 /* Register the tile */ 640 #pragma omp critical 641 list_add_tail(&tiles, &tile->node); 642 643 /* Setup the tile coordinates */ 644 tile->data.x = (uint16_t)tile_org[0]; 645 tile->data.y = (uint16_t)tile_org[1]; 646 647 /* Setup the tile coordinates in the image plane */ 648 tile_org[0] *= TILE_SIZE; 649 tile_org[1] *= TILE_SIZE; 650 tile_sz[0] = MMIN(TILE_SIZE, args->image_definition[0] - tile_org[0]); 651 tile_sz[1] = MMIN(TILE_SIZE, args->image_definition[1] - tile_org[1]); 652 653 /* Draw the tile */ 654 res_local = solve_tile 655 (scn, rng, enc_id, args->cam, args->time_range, tile_org, tile_sz, 656 args->spp, register_paths, pix_sz, args->picard_order, args->diff_algo, 657 buf, tile); 658 if(res_local != RES_OK) { 659 ATOMIC_SET(&res, res_local); 660 continue; 661 } 662 663 /* Update progress */ 664 n = (size_t)ATOMIC_INCR(&nsolved_tiles); 665 pcent = (int)((double)n*100.0 / (double)ntiles_proc + 0.5/*round*/); 666 #pragma omp critical 667 if(pcent/pcent_progress > progress[0]/pcent_progress) { 668 progress[0] = pcent; 669 print_progress_update(scn->dev, progress, PROGRESS_MSG); 670 } 671 } 672 673 /* Synchronise the processes */ 674 process_barrier(scn->dev); 675 676 res = gather_res_T(scn->dev, (res_T)res); 677 if(res != RES_OK) goto error; 678 679 print_progress_completion(scn->dev, progress, PROGRESS_MSG); 680 #undef PROGRESS_MSG 681 682 /* Report computation time */ 683 time_sub(&time0, time_current(&time1), &time0); 684 time_dump(&time0, TIME_ALL, NULL, buffer, sizeof(buffer)); 685 log_info(scn->dev, "Image rendered in %s.\n", buffer); 686 687 /* Gather the RNG proxy sequence IDs and ensure that the RNG proxy state of 688 * the master process is greater than the RNG proxy state of all other 689 * processes */ 690 res = gather_rng_proxy_sequence_id(scn->dev, rng_proxy); 691 if(res != RES_OK) goto error; 692 693 time_current(&time0); 694 695 res = gather_tiles(scn->dev, buf, args->spp, ntiles, &tiles); 696 if(res != RES_OK) goto error; 697 698 time_sub(&time0, time_current(&time1), &time0); 699 time_dump(&time0, TIME_ALL, NULL, buffer, sizeof(buffer)); 700 log_info(scn->dev, "Image tiles gathered in %s.\n", buffer); 701 702 if(is_master_process) { 703 res = finalize_estimator_buffer(buf, rng_proxy, args->spp); 704 if(res != RES_OK) goto error; 705 } 706 707 exit: 708 free_rendered_tiles(&tiles); 709 if(per_thread_rng)release_per_thread_rng(scn->dev, per_thread_rng); 710 if(progress) free_process_progress(scn->dev, progress); 711 if(rng_proxy) SSP(rng_proxy_ref_put(rng_proxy)); 712 if(out_buf) *out_buf = buf; 713 return (res_T)res; 714 error: 715 if(buf) { SDIS(estimator_buffer_ref_put(buf)); buf = NULL; } 716 goto exit; 717 }