htrdr_planets.c (24333B)
1 /* Copyright (C) 2018-2019, 2022-2025 Centre National de la Recherche Scientifique 2 * Copyright (C) 2020-2022 Institut Mines Télécom Albi-Carmaux 3 * Copyright (C) 2022-2025 Institut Pierre-Simon Laplace 4 * Copyright (C) 2022-2025 Institut de Physique du Globe de Paris 5 * Copyright (C) 2018-2025 |Méso|Star> (contact@meso-star.com) 6 * Copyright (C) 2022-2025 Observatoire de Paris 7 * Copyright (C) 2022-2025 Université de Reims Champagne-Ardenne 8 * Copyright (C) 2022-2025 Université de Versaille Saint-Quentin 9 * Copyright (C) 2018-2019, 2022-2025 Université Paul Sabatier 10 * 11 * This program is free software: you can redistribute it and/or modify 12 * it under the terms of the GNU General Public License as published by 13 * the Free Software Foundation, either version 3 of the License, or 14 * (at your option) any later version. 15 * 16 * This program is distributed in the hope that it will be useful, 17 * but WITHOUT ANY WARRANTY; without even the implied warranty of 18 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 19 * GNU General Public License for more details. 20 * 21 * You should have received a copy of the GNU General Public License 22 * along with this program. If not, see <http://www.gnu.org/licenses/>. */ 23 24 #define _POSIX_C_SOURCE 200112L /* fdopen, nanosleep, nextafter, rint */ 25 26 #include "core/htrdr.h" 27 #include "core/htrdr_ran_wlen_cie_xyz.h" 28 #include "core/htrdr_ran_wlen_discrete.h" 29 #include "core/htrdr_ran_wlen_planck.h" 30 #include "core/htrdr_log.h" 31 32 #include "planets/htrdr_planets.h" 33 #include "planets/htrdr_planets_args.h" 34 #include "planets/htrdr_planets_c.h" 35 #include "planets/htrdr_planets_source.h" 36 37 #include <rad-net/rnatm.h> 38 #include <rad-net/rngrd.h> 39 40 #include <star/scam.h> 41 #include <star/smsh.h> 42 43 #include <rsys/cstr.h> 44 #include <rsys/double3.h> 45 #include <rsys/mem_allocator.h> 46 47 #include <fcntl.h> /* open */ 48 #include <math.h> /* nextafter, rint */ 49 #include <unistd.h> /* close */ 50 #include <sys/stat.h> 51 52 #include <mpi.h> 53 #include <time.h> 54 55 /******************************************************************************* 56 * Helper function 57 ******************************************************************************/ 58 /* Calculate the number of fixed size spectral intervals to use for the 59 * cumulative */ 60 static size_t 61 compute_nintervals_for_spectral_cdf(const struct htrdr_planets* cmd) 62 { 63 double range_size; 64 size_t nintervals; 65 ASSERT(cmd); 66 67 range_size = 68 cmd->spectral_domain.wlen_range[1] 69 - cmd->spectral_domain.wlen_range[0]; 70 71 /* Initially assume ~one interval per nanometer */ 72 nintervals = (size_t)rint(range_size); 73 74 return nintervals; 75 } 76 77 static res_T 78 setup_octree_storage 79 (struct htrdr_planets* cmd, 80 const struct htrdr_planets_args* args, 81 struct rnatm_create_args* rnatm_args) 82 { 83 struct stat file_stat; 84 int fd = -1; 85 int err = 0; 86 res_T res = RES_OK; 87 ASSERT(cmd && args && rnatm_args); 88 89 rnatm_args->octrees_storage = NULL; 90 rnatm_args->load_octrees_from_storage = 0; 91 92 if(!args->accel_struct.storage) goto exit; 93 94 fd = open(args->accel_struct.storage, O_CREAT|O_RDWR, S_IRUSR|S_IWUSR); 95 if(fd < 0) { res = RES_IO_ERR; goto error; } 96 97 rnatm_args->octrees_storage = fdopen(fd, "w+"); 98 if(!rnatm_args->octrees_storage) { res = RES_IO_ERR; goto error; } 99 100 /* From now on, manage the opened file from its pointer and not from its 101 * descriptor */ 102 fd = -1; 103 104 err = stat(args->accel_struct.storage, &file_stat); 105 if(err < 0) { res = RES_IO_ERR; goto error; } 106 107 if(file_stat.st_size != 0) { 108 /* The file is not empty and therefore must contain valid octrees */ 109 rnatm_args->load_octrees_from_storage = 1; 110 } 111 112 exit: 113 cmd->octrees_storage = rnatm_args->octrees_storage; 114 return res; 115 error: 116 htrdr_log_err(cmd->htrdr, "error opening the octree storage `%s' -- %s\n", 117 args->accel_struct.storage, res_to_cstr(res)); 118 119 if(fd >= 0) CHK(close(fd) == 0); 120 if(rnatm_args->octrees_storage) CHK(fclose(rnatm_args->octrees_storage) == 0); 121 rnatm_args->octrees_storage = NULL; 122 rnatm_args->load_octrees_from_storage = 1; 123 goto exit; 124 } 125 126 static void 127 mpi_barrier(void) 128 { 129 struct timespec t; 130 int complete = 0; 131 MPI_Request req; 132 133 t.tv_sec = 0; 134 t.tv_nsec = 10000000; /* 10ms */ 135 136 /* Use an asynchronous barrier to avoid active waiting, 137 * and therefore wasted resources */ 138 MPI_Ibarrier(MPI_COMM_WORLD, &req); 139 140 do { 141 nanosleep(&t, NULL); 142 MPI_Test(&req, &complete, MPI_STATUS_IGNORE); 143 } while(!complete); 144 } 145 146 static res_T 147 setup_atmosphere 148 (struct htrdr_planets* cmd, 149 const struct htrdr_planets_args* args) 150 { 151 struct rnatm_create_args rnatm = RNATM_CREATE_ARGS_DEFAULT; 152 res_T res = RES_OK; 153 ASSERT(cmd && args); 154 155 rnatm.gas = args->gas; 156 rnatm.aerosols = args->aerosols; 157 rnatm.naerosols = args->naerosols; 158 rnatm.name = "atmosphere"; 159 rnatm.spectral_range[0] = args->spectral_domain.wlen_range[0]; 160 rnatm.spectral_range[1] = args->spectral_domain.wlen_range[1]; 161 rnatm.optical_thickness = args->accel_struct.optical_thickness; 162 rnatm.grid_definition_hint = args->accel_struct.definition_hint; 163 rnatm.precompute_normals = args->precompute_normals; 164 rnatm.logger = htrdr_get_logger(cmd->htrdr); 165 rnatm.allocator = htrdr_get_allocator(cmd->htrdr); 166 rnatm.nthreads = args->accel_struct.nthreads; 167 rnatm.verbose = args->verbose; 168 169 if(!args->accel_struct.storage || !args->accel_struct.master_only) { 170 /* From now on, either the octrees have to be built in memory, or all the 171 * processes have to build them. In all cases, all processes must create 172 * the atmopshere data structures. */ 173 if((res = setup_octree_storage(cmd, args, &rnatm)) != RES_OK) goto error; 174 if((res = rnatm_create(&rnatm, &cmd->atmosphere)) != RES_OK) goto error; 175 176 } else { 177 178 const int is_master_process = htrdr_get_mpi_rank(cmd->htrdr) == 0; 179 /* A storage space is used and only the master process must fill it with the 180 * octrees it builds */ 181 182 if(is_master_process) { 183 /* Octrees are built only by the master process and stored on disk. 184 * Note that in reality, octrees may already be constructed and stored in 185 * the storage provided. In any case, we can treat this special case as 186 * the general case, and therefore simply consider that in such situation, 187 * "construction" by the master process is only faster */ 188 if((res = setup_octree_storage(cmd, args, &rnatm)) != RES_OK) goto error; 189 if((res = rnatm_create(&rnatm, &cmd->atmosphere)) != RES_OK) goto error; 190 } 191 192 /* Wait for octrees to be built by the master process */ 193 mpi_barrier(); 194 195 if(!is_master_process) { 196 /* Octrees have been built by the master process. 197 * Non master processes simply load them. */ 198 if((res = setup_octree_storage(cmd, args, &rnatm)) != RES_OK) goto error; 199 ASSERT(rnatm.load_octrees_from_storage == 1); 200 if((res = rnatm_create(&rnatm, &cmd->atmosphere)) != RES_OK) goto error; 201 } 202 } 203 204 exit: 205 return res; 206 error: 207 if(cmd->atmosphere) { 208 RNATM(ref_put(cmd->atmosphere)); 209 cmd->atmosphere = NULL; 210 } 211 goto exit; 212 } 213 214 static res_T 215 setup_ground 216 (struct htrdr_planets* cmd, 217 const struct htrdr_planets_args* args) 218 { 219 struct rngrd_create_args rngrd_args = RNGRD_CREATE_ARGS_DEFAULT; 220 res_T res = RES_OK; 221 ASSERT(cmd && args); 222 223 if(cmd->output_type == HTRDR_PLANETS_ARGS_OUTPUT_OCTREES) 224 goto exit; 225 226 rngrd_args.smsh_filename = args->ground.smsh_filename; 227 rngrd_args.props_filename = args->ground.props_filename; 228 rngrd_args.mtllst_filename = args->ground.mtllst_filename; 229 rngrd_args.name = args->ground.name; 230 rngrd_args.logger = htrdr_get_logger(cmd->htrdr); 231 rngrd_args.allocator = htrdr_get_allocator(cmd->htrdr); 232 rngrd_args.verbose = args->verbose; 233 234 res = rngrd_create(&rngrd_args, &cmd->ground); 235 if(res != RES_OK) goto error; 236 237 exit: 238 return res; 239 error: 240 if(cmd->ground) { 241 RNGRD(ref_put(cmd->ground)); 242 cmd->ground = NULL; 243 } 244 goto exit; 245 } 246 247 static res_T 248 setup_spectral_domain_sw 249 (struct htrdr_planets* cmd, 250 const struct htrdr_planets_args* args) 251 { 252 res_T res = RES_OK; 253 ASSERT(cmd && args); 254 ASSERT(cmd->spectral_domain.type == HTRDR_SPECTRAL_SW); 255 256 /* Discrete distribution */ 257 if(args->source.rnrl_filename) { 258 struct htrdr_planets_source_spectrum spectrum; 259 struct htrdr_ran_wlen_discrete_create_args discrete_args; 260 261 res = htrdr_planets_source_get_spectrum 262 (cmd->source, cmd->spectral_domain.wlen_range, &spectrum); 263 if(res != RES_OK) goto error; 264 265 discrete_args.get = htrdr_planets_source_spectrum_at; 266 discrete_args.nwavelengths = spectrum.size; 267 discrete_args.context = &spectrum; 268 res = htrdr_ran_wlen_discrete_create 269 (cmd->htrdr, &discrete_args, &cmd->discrete); 270 if(res != RES_OK) goto error; 271 272 /* Planck distribution */ 273 } else { 274 const size_t nintervals = compute_nintervals_for_spectral_cdf(cmd); 275 276 /* Use the source temperature as the reference temperature of the Planck 277 * distribution */ 278 res = htrdr_ran_wlen_planck_create(cmd->htrdr, 279 cmd->spectral_domain.wlen_range, nintervals, args->source.temperature, 280 &cmd->planck); 281 if(res != RES_OK) goto error; 282 } 283 284 exit: 285 return res; 286 error: 287 goto exit; 288 } 289 290 static INLINE res_T 291 setup_spectral_domain 292 (struct htrdr_planets* cmd, 293 const struct htrdr_planets_args* args) 294 { 295 double ground_T_range[2]; 296 size_t nintervals; 297 res_T res = RES_OK; 298 ASSERT(cmd && args); 299 300 /* No spectral distribution required to write octrees */ 301 if(cmd->output_type == HTRDR_PLANETS_ARGS_OUTPUT_OCTREES) 302 goto exit; 303 304 cmd->spectral_domain = args->spectral_domain; 305 306 /* Configure the spectral distribution */ 307 switch(cmd->spectral_domain.type) { 308 309 case HTRDR_SPECTRAL_LW: 310 res = rngrd_get_temperature_range(cmd->ground, ground_T_range); 311 if(res != RES_OK) goto error; 312 313 /* Use as the reference temperature of the Planck distribution the 314 * maximum scene temperature which, in fact, should be the maximum ground 315 * temperature */ 316 nintervals = compute_nintervals_for_spectral_cdf(cmd); 317 res = htrdr_ran_wlen_planck_create(cmd->htrdr, 318 cmd->spectral_domain.wlen_range, nintervals, ground_T_range[1], 319 &cmd->planck); 320 if(res != RES_OK) goto error; 321 break; 322 323 case HTRDR_SPECTRAL_SW: 324 res = setup_spectral_domain_sw(cmd, args); 325 if(res != RES_OK) goto error; 326 break; 327 328 case HTRDR_SPECTRAL_SW_CIE_XYZ: 329 /* CIE XYZ distribution */ 330 nintervals = compute_nintervals_for_spectral_cdf(cmd); 331 res = htrdr_ran_wlen_cie_xyz_create(cmd->htrdr, 332 cmd->spectral_domain.wlen_range, nintervals, &cmd->cie); 333 if(res != RES_OK) goto error; 334 break; 335 336 default: FATAL("Unreachable code\n"); break; 337 } 338 339 exit: 340 return res; 341 error: 342 goto exit; 343 } 344 345 static res_T 346 setup_output 347 (struct htrdr_planets* cmd, 348 const struct htrdr_planets_args* args) 349 { 350 const char* output_name = NULL; 351 res_T res = RES_OK; 352 ASSERT(cmd && args); 353 354 /* No output stream on non master processes */ 355 if(htrdr_get_mpi_rank(cmd->htrdr) != 0) { 356 cmd->output = NULL; 357 output_name = "<null>"; 358 359 /* Write results on stdout */ 360 } else if(!args->output) { 361 cmd->output = stdout; 362 output_name = "<stdout>"; 363 364 /* Open the output stream */ 365 } else { 366 res = htrdr_open_output_stream(cmd->htrdr, args->output, 0/*read*/, 367 args->force_output_overwrite, &cmd->output); 368 if(res != RES_OK) goto error; 369 output_name = args->output; 370 } 371 372 res = str_set(&cmd->output_name, output_name); 373 if(res != RES_OK) { 374 htrdr_log_err(cmd->htrdr, "error storing output stream name `%s' -- %s\n", 375 output_name, res_to_cstr(res)); 376 goto error; 377 } 378 379 cmd->output_type = args->output_type; 380 381 exit: 382 return res; 383 error: 384 str_clear(&cmd->output_name); 385 if(cmd->output && cmd->output != stdout) { 386 CHK(fclose(cmd->output) == 0); 387 cmd->output = NULL; 388 } 389 goto exit; 390 } 391 392 static INLINE res_T 393 setup_source 394 (struct htrdr_planets* cmd, 395 const struct htrdr_planets_args* args) 396 { 397 res_T res = RES_OK; 398 ASSERT(cmd && args); 399 400 if(cmd->output_type == HTRDR_PLANETS_ARGS_OUTPUT_OCTREES 401 /* No source in Longwave. 402 * Check the spectral domain type on the input arguments, because when the 403 * source is configured, the spectral dimension may not yet be set */ 404 || args->spectral_domain.type == HTRDR_SPECTRAL_LW) 405 goto exit; 406 407 res = htrdr_planets_source_create(cmd->htrdr, &args->source, &cmd->source); 408 if(res != RES_OK) goto error; 409 410 exit: 411 return res; 412 error: 413 goto exit; 414 } 415 416 static res_T 417 setup_camera_orthographic 418 (struct htrdr_planets* cmd, 419 const struct htrdr_planets_args* args) 420 { 421 struct scam_orthographic_args cam_args = SCAM_ORTHOGRAPHIC_ARGS_DEFAULT; 422 res_T res = RES_OK; 423 424 /* Check pre-conditions */ 425 ASSERT(cmd && args); 426 ASSERT(cmd->output_type == HTRDR_PLANETS_ARGS_OUTPUT_IMAGE); 427 ASSERT(args->cam_type == HTRDR_ARGS_CAMERA_ORTHOGRAPHIC); 428 ASSERT(htrdr_args_camera_orthographic_check(&args->cam_ortho) == RES_OK); 429 ASSERT(htrdr_args_image_check(&args->image) == RES_OK); 430 431 d3_set(cam_args.position, args->cam_ortho.position); 432 d3_set(cam_args.target, args->cam_ortho.target); 433 d3_set(cam_args.up, args->cam_ortho.up); 434 cam_args.height = args->cam_ortho.height; 435 cam_args.aspect_ratio = 436 (double)args->image.definition[0] 437 / (double)args->image.definition[1]; 438 439 res = scam_create_orthographic 440 (htrdr_get_logger(cmd->htrdr), 441 htrdr_get_allocator(cmd->htrdr), 442 htrdr_get_verbosity_level(cmd->htrdr), 443 &cam_args, 444 &cmd->camera); 445 if(res != RES_OK) goto error; 446 447 exit: 448 return res; 449 error: 450 goto exit; 451 } 452 453 static res_T 454 setup_camera_perspective 455 (struct htrdr_planets* cmd, 456 const struct htrdr_planets_args* args) 457 { 458 struct scam_perspective_args cam_args = SCAM_PERSPECTIVE_ARGS_DEFAULT; 459 res_T res = RES_OK; 460 461 /* Check pre-conditions */ 462 ASSERT(cmd && args); 463 ASSERT(args->output_type == HTRDR_PLANETS_ARGS_OUTPUT_IMAGE); 464 ASSERT(args->cam_type == HTRDR_ARGS_CAMERA_PERSPECTIVE); 465 ASSERT(htrdr_args_camera_perspective_check(&args->cam_persp) == RES_OK); 466 ASSERT(htrdr_args_image_check(&args->image) == RES_OK); 467 468 d3_set(cam_args.position, args->cam_persp.position); 469 d3_set(cam_args.target, args->cam_persp.target); 470 d3_set(cam_args.up, args->cam_persp.up); 471 cam_args.field_of_view = MDEG2RAD(args->cam_persp.fov_y); 472 cam_args.lens_radius = args->cam_persp.lens_radius; 473 cam_args.focal_distance = args->cam_persp.focal_dst; 474 cam_args.aspect_ratio = 475 (double)args->image.definition[0] 476 / (double)args->image.definition[1]; 477 478 res = scam_create_perspective 479 (htrdr_get_logger(cmd->htrdr), 480 htrdr_get_allocator(cmd->htrdr), 481 htrdr_get_verbosity_level(cmd->htrdr), 482 &cam_args, 483 &cmd->camera); 484 if(res != RES_OK) goto error; 485 486 exit: 487 return res; 488 error: 489 goto exit; 490 } 491 492 static res_T 493 setup_camera 494 (struct htrdr_planets* cmd, 495 const struct htrdr_planets_args* args) 496 { 497 res_T res = RES_OK; 498 ASSERT(cmd->output_type == HTRDR_PLANETS_ARGS_OUTPUT_IMAGE); 499 500 switch(args->cam_type) { 501 case HTRDR_ARGS_CAMERA_ORTHOGRAPHIC: 502 res = setup_camera_orthographic(cmd, args); 503 break; 504 case HTRDR_ARGS_CAMERA_PERSPECTIVE: 505 res = setup_camera_perspective(cmd, args); 506 break; 507 default: FATAL("Unreachable code.\n"); break; 508 } 509 510 return res; 511 } 512 513 static res_T 514 setup_buffer_image 515 (struct htrdr_planets* cmd, 516 const struct htrdr_planets_args* args) 517 { 518 struct htrdr_pixel_format pixfmt = HTRDR_PIXEL_FORMAT_NULL; 519 res_T res = RES_OK; 520 ASSERT(cmd && args); 521 ASSERT(cmd->output_type == HTRDR_PLANETS_ARGS_OUTPUT_IMAGE); 522 523 planets_get_pixel_format(cmd, &pixfmt); 524 525 /* Setup buffer layout */ 526 cmd->buf_layout.width = args->image.definition[0]; 527 cmd->buf_layout.height = args->image.definition[1]; 528 cmd->buf_layout.pitch = args->image.definition[0] * pixfmt.size; 529 cmd->buf_layout.elmt_size = pixfmt.size; 530 cmd->buf_layout.alignment = pixfmt.alignment; 531 532 /* Save the number of samples per pixel */ 533 cmd->spp = args->image.spp; 534 535 /* Create the image buffer only on the master process; Image parts rendered 536 * by other processes are collected there */ 537 if(htrdr_get_mpi_rank(cmd->htrdr) != 0) goto exit; 538 539 res = htrdr_buffer_create(cmd->htrdr, &cmd->buf_layout, &cmd->buf); 540 if(res != RES_OK) goto error; 541 542 exit: 543 return res; 544 error: 545 if(cmd->buf) { htrdr_buffer_ref_put(cmd->buf); cmd->buf = NULL; } 546 goto exit; 547 } 548 549 static res_T 550 setup_buffer_raw 551 (struct htrdr_planets* cmd, 552 const struct htrdr_planets_args* args) 553 { 554 struct smsh_desc desc = SMSH_DESC_NULL; 555 size_t sz = 0; /* Size of a voxel storing volumic radiative budget */ 556 size_t al = 0; /* Alignment of a voxel storing volumic radiative budget */ 557 res_T res = RES_OK; 558 559 ASSERT(cmd && args); 560 ASSERT(cmd->output_type == HTRDR_PLANETS_ARGS_OUTPUT_VOLUMIC_RADIATIVE_BUDGET); 561 ASSERT(cmd->volrad_mesh != NULL); /* The volurad mesh must be defined */ 562 563 res = smsh_get_desc(cmd->volrad_mesh, &desc); 564 if(res != RES_OK) goto error; 565 566 /* Setup buffer layout for volumic radiative budget calculation */ 567 sz = sizeof(struct planets_voxel_radiative_budget); 568 al = ALIGNOF(struct planets_voxel_radiative_budget); 569 cmd->buf_layout.width = desc.ncells; 570 cmd->buf_layout.height = 1; 571 cmd->buf_layout.pitch = desc.ncells * sz; 572 cmd->buf_layout.elmt_size = sz; 573 cmd->buf_layout.alignment = al; 574 575 /* Save the number of samples per tetrahedron */ 576 cmd->spt = args->volrad_budget.spt; 577 578 /* Create the raw buffer only on master process; buffer parts calculated by 579 * other processes are collected there */ 580 if(htrdr_get_mpi_rank(cmd->htrdr) != 0) goto exit; 581 582 res = htrdr_buffer_create(cmd->htrdr, &cmd->buf_layout, &cmd->buf); 583 if(res != RES_OK) goto error; 584 585 exit: 586 return res; 587 error: 588 if(cmd->buf) { htrdr_buffer_ref_put(cmd->buf); cmd->buf = NULL; } 589 goto exit; 590 } 591 592 static res_T 593 setup_buffer 594 (struct htrdr_planets* cmd, 595 const struct htrdr_planets_args* args) 596 { 597 res_T res = RES_OK; 598 ASSERT(cmd && args); 599 600 switch(cmd->output_type) { 601 case HTRDR_PLANETS_ARGS_OUTPUT_IMAGE: 602 res = setup_buffer_image(cmd, args); 603 break; 604 case HTRDR_PLANETS_ARGS_OUTPUT_VOLUMIC_RADIATIVE_BUDGET: 605 res = setup_buffer_raw(cmd, args); 606 break; 607 default: /* Nothing to do */ break; 608 } 609 if(res != RES_OK) goto error; 610 611 exit: 612 return res; 613 error: 614 goto exit; 615 } 616 617 static res_T 618 setup_volrad_budget_mesh 619 (struct htrdr_planets* cmd, 620 const struct htrdr_planets_args* args) 621 { 622 struct smsh_create_args create_args = SMSH_CREATE_ARGS_DEFAULT; 623 struct smsh_load_args load_args = SMSH_LOAD_ARGS_NULL; 624 struct smsh_desc desc = SMSH_DESC_NULL; 625 res_T res = RES_OK; 626 ASSERT(cmd && args); 627 628 if(cmd->output_type != HTRDR_PLANETS_ARGS_OUTPUT_VOLUMIC_RADIATIVE_BUDGET) 629 goto exit; 630 631 /* Store the number of samples per tetrahedron to be used */ 632 cmd->spt = args->volrad_budget.spt; 633 634 create_args.logger = htrdr_get_logger(cmd->htrdr); 635 create_args.allocator = htrdr_get_allocator(cmd->htrdr); 636 create_args.verbose = htrdr_get_verbosity_level(cmd->htrdr); 637 res = smsh_create(&create_args, &cmd->volrad_mesh); 638 if(res != RES_OK) goto error; 639 640 load_args.path = args->volrad_budget.smsh_filename; 641 res = smsh_load(cmd->volrad_mesh, &load_args); 642 if(res != RES_OK) goto error; 643 644 res = smsh_get_desc(cmd->volrad_mesh, &desc); 645 if(res != RES_OK) goto error; 646 647 /* Check that the loaded mesh is effectively a volume mesh */ 648 if(desc.dnode != 3 || desc.dcell != 4) { 649 htrdr_log_err(cmd->htrdr, 650 "%s: the volumic radiative budget calculation " 651 "expects a 3D tetrahedral mesh " 652 "(dimension of mesh: %u; dimension of the vertices: %u)\n", 653 args->volrad_budget.smsh_filename, 654 desc.dnode, 655 desc.dcell); 656 res = RES_BAD_ARG; 657 goto error; 658 } 659 660 exit: 661 return res; 662 error: 663 if(cmd->volrad_mesh) { 664 SMSH(ref_put(cmd->volrad_mesh)); 665 cmd->volrad_mesh = NULL; 666 } 667 goto exit; 668 } 669 670 static INLINE res_T 671 write_vtk_octrees(const struct htrdr_planets* cmd) 672 { 673 size_t octrees_range[2]; 674 res_T res = RES_OK; 675 ASSERT(cmd); 676 677 /* Nothing to do on non master process */ 678 if(htrdr_get_mpi_rank(cmd->htrdr) != 0) goto exit; 679 680 octrees_range[0] = 0; 681 octrees_range[1] = rnatm_get_spectral_items_count(cmd->atmosphere) - 1; 682 683 res = rnatm_write_vtk_octrees(cmd->atmosphere, octrees_range, cmd->output); 684 if(res != RES_OK) goto error; 685 686 exit: 687 return res; 688 error: 689 goto exit; 690 } 691 692 static void 693 planets_release(ref_T* ref) 694 { 695 struct htrdr_planets* cmd = CONTAINER_OF(ref, struct htrdr_planets, ref); 696 struct htrdr* htrdr = NULL; 697 ASSERT(ref); 698 699 if(cmd->atmosphere) RNATM(ref_put(cmd->atmosphere)); 700 if(cmd->ground) RNGRD(ref_put(cmd->ground)); 701 if(cmd->source) htrdr_planets_source_ref_put(cmd->source); 702 if(cmd->cie) htrdr_ran_wlen_cie_xyz_ref_put(cmd->cie); 703 if(cmd->discrete) htrdr_ran_wlen_discrete_ref_put(cmd->discrete); 704 if(cmd->planck) htrdr_ran_wlen_planck_ref_put(cmd->planck); 705 if(cmd->octrees_storage) CHK(fclose(cmd->octrees_storage) == 0); 706 if(cmd->output && cmd->output != stdout) CHK(fclose(cmd->output) == 0); 707 if(cmd->buf) htrdr_buffer_ref_put(cmd->buf); 708 if(cmd->camera) SCAM(ref_put(cmd->camera)); 709 if(cmd->volrad_mesh) SMSH(ref_put(cmd->volrad_mesh)); 710 str_release(&cmd->output_name); 711 712 htrdr = cmd->htrdr; 713 MEM_RM(htrdr_get_allocator(htrdr), cmd); 714 htrdr_ref_put(htrdr); 715 } 716 717 /******************************************************************************* 718 * Exported functions 719 ******************************************************************************/ 720 res_T 721 htrdr_planets_create 722 (struct htrdr* htrdr, 723 const struct htrdr_planets_args* args, 724 struct htrdr_planets** out_cmd) 725 { 726 struct htrdr_planets* cmd = NULL; 727 res_T res = RES_OK; 728 ASSERT(htrdr && out_cmd); 729 730 res = htrdr_planets_args_check(args); 731 if(res != RES_OK) { 732 htrdr_log_err(htrdr, "Invalid htrdr_planets arguments -- %s\n", 733 res_to_cstr(res)); 734 goto error; 735 } 736 737 cmd = MEM_CALLOC(htrdr_get_allocator(htrdr), 1, sizeof(*cmd)); 738 if(!cmd) { 739 htrdr_log_err(htrdr, "Error allocating htrdr_planets command\n"); 740 res = RES_MEM_ERR; 741 goto error; 742 } 743 ref_init(&cmd->ref); 744 htrdr_ref_get(htrdr); 745 cmd->htrdr = htrdr; 746 str_init(htrdr_get_allocator(htrdr), &cmd->output_name); 747 748 res = setup_output(cmd, args); 749 if(res != RES_OK) goto error; 750 res = setup_source(cmd, args); 751 if(res != RES_OK) goto error; 752 res = setup_camera(cmd, args); 753 if(res != RES_OK) goto error; 754 res = setup_ground(cmd, args); 755 if(res != RES_OK) goto error; 756 res = setup_atmosphere(cmd, args); 757 if(res != RES_OK) goto error; 758 res = setup_spectral_domain(cmd, args); 759 if(res != RES_OK) goto error; 760 res = setup_volrad_budget_mesh(cmd, args); 761 if(res != RES_OK) goto error; 762 res = setup_buffer(cmd, args); 763 if(res != RES_OK) goto error; 764 765 exit: 766 *out_cmd = cmd; 767 return res; 768 error: 769 if(cmd) { 770 htrdr_planets_ref_put(cmd); 771 cmd = NULL; 772 } 773 goto exit; 774 } 775 776 void 777 htrdr_planets_ref_get(struct htrdr_planets* cmd) 778 { 779 ASSERT(cmd); 780 ref_get(&cmd->ref); 781 } 782 783 void 784 htrdr_planets_ref_put(struct htrdr_planets* cmd) 785 { 786 ASSERT(cmd); 787 ref_put(&cmd->ref, planets_release); 788 } 789 790 res_T 791 htrdr_planets_run(struct htrdr_planets* cmd) 792 { 793 res_T res = RES_OK; 794 ASSERT(cmd); 795 796 switch(cmd->output_type) { 797 case HTRDR_PLANETS_ARGS_OUTPUT_IMAGE: 798 res = planets_draw_map(cmd); 799 break; 800 case HTRDR_PLANETS_ARGS_OUTPUT_OCTREES: 801 res = write_vtk_octrees(cmd); 802 break; 803 case HTRDR_PLANETS_ARGS_OUTPUT_VOLUMIC_RADIATIVE_BUDGET: 804 res = planets_solve_volrad_budget(cmd); 805 break; 806 default: FATAL("Unreachable code\n"); break; 807 } 808 if(res != RES_OK) goto error; 809 810 exit: 811 return res; 812 error: 813 goto exit; 814 } 815 816 /******************************************************************************* 817 * Local function 818 ******************************************************************************/ 819 void 820 planets_get_pixel_format 821 (const struct htrdr_planets* cmd, 822 struct htrdr_pixel_format* fmt) 823 { 824 ASSERT(cmd && fmt && cmd->output_type == HTRDR_PLANETS_ARGS_OUTPUT_IMAGE); 825 (void)cmd; 826 827 switch(cmd->spectral_domain.type) { 828 case HTRDR_SPECTRAL_LW: 829 case HTRDR_SPECTRAL_SW: 830 fmt->size = sizeof(struct planets_pixel_xwave); 831 fmt->alignment = ALIGNOF(struct planets_pixel_xwave); 832 break; 833 case HTRDR_SPECTRAL_SW_CIE_XYZ: 834 fmt->size = sizeof(struct planets_pixel_image); 835 fmt->alignment = ALIGNOF(struct planets_pixel_image); 836 break; 837 default: FATAL("Unreachable code\n"); break; 838 } 839 }