sln_get.c (14427B)
1 /* Copyright (C) 2022, 2026 |Méso|Star> (contact@meso-star.com) 2 * Copyright (C) 2026 Université de Lorraine 3 * Copyright (C) 2022 Centre National de la Recherche Scientifique 4 * Copyright (C) 2022 Université Paul Sabatier 5 * 6 * This program is free software: you can redistribute it and/or modify 7 * it under the terms of the GNU General Public License as published by 8 * the Free Software Foundation, either version 3 of the License, or 9 * (at your option) any later version. 10 * 11 * This program is distributed in the hope that it will be useful, 12 * but WITHOUT ANY WARRANTY; without even the implied warranty of 13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 14 * GNU General Public License for more details. 15 * 16 * You should have received a copy of the GNU General Public License 17 * along with this program. If not, see <http://www.gnu.org/licenses/>. */ 18 19 #define _POSIX_C_SOURCE 200112L /* getopt */ 20 21 #include "sln.h" 22 23 #include <rsys/cstr.h> 24 #include <rsys/math.h> 25 #include <rsys/mem_allocator.h> 26 27 #include <unistd.h> 28 29 enum child { LEFT, RIGHT }; 30 31 enum output_type { 32 OUTPUT_LEVEL_DESCRIPTOR, 33 OUTPUT_NODE_DESCRIPTOR, 34 OUTPUT_NODE_MESH, 35 OUTPUT_NODE_VALUE, 36 OUTPUT_TREE_DESCRIPTOR, 37 OUTPUT_COUNT__ 38 }; 39 40 struct args { 41 const char* tree; /* NULL <=> read from standard input */ 42 const char* molparams; 43 const char* lines; 44 45 enum output_type output_type; 46 47 double wavenumber; /* Wave number at which the spectrum is evaluated */ 48 49 /* Steps for traversing the tree */ 50 unsigned descent_path[SLN_TREE_DEPTH_MAX]; 51 unsigned depth; /* Current depth in the tree */ 52 53 /* Miscellaneous */ 54 unsigned level; /* Queried level */ 55 int quit; 56 int verbose; 57 int lines_in_shtr_format; 58 }; 59 #define ARGS_DEFAULT__ {NULL,NULL,NULL,OUTPUT_TREE_DESCRIPTOR,0,{0},0,0,0,0,0} 60 static const struct args ARGS_DEFAULT = ARGS_DEFAULT__; 61 62 struct cmd { 63 struct args args; 64 65 struct sln_device* sln; 66 struct sln_tree* tree; 67 68 struct shtr* shtr; 69 struct shtr_isotope_metadata* molparams; 70 struct shtr_line_list* lines; 71 }; 72 #define CMD_NULL__ {0} 73 static const struct cmd CMD_NULL = CMD_NULL__; 74 75 /******************************************************************************* 76 * Helper functions 77 ******************************************************************************/ 78 static void 79 usage(FILE* stream) 80 { 81 fprintf(stream, 82 "usage: sln-get [-hlmnrsv] [-c child_id[:level_count]] [-d level] [-w wavenumber]\n" 83 " -i lines -p molparams [tree]\n"); 84 } 85 86 static res_T 87 tree_descent(struct args* args, const char* str) 88 { 89 unsigned path[2] = {0/*child_id*/,1/*level_count*/}; 90 unsigned i=0; 91 size_t n=0; 92 res_T res = RES_OK; 93 ASSERT(args && str); 94 95 res = cstr_to_list_uint(str, ':', path, &n, 2); 96 if(res != RES_OK) goto error; 97 98 for(i=0; i<path[1] && args->depth<SLN_TREE_DEPTH_MAX; ++i, ++args->depth) { 99 args->descent_path[args->depth] = path[0]; 100 } 101 102 exit: 103 return res; 104 error: 105 goto exit; 106 } 107 108 static res_T 109 args_init(struct args* args, int argc, char** argv) 110 { 111 int opt = 0; 112 res_T res = RES_OK; 113 114 ASSERT(args); 115 116 *args = ARGS_DEFAULT; 117 118 while((opt = getopt(argc, argv, "c:d:hi:mnp:svw:")) != -1) { 119 switch(opt) { 120 case 'c': res = tree_descent(args, optarg); break; 121 case 'd': 122 args->output_type = OUTPUT_LEVEL_DESCRIPTOR; 123 res = cstr_to_uint(optarg, &args->level); 124 break; 125 case 'h': 126 usage(stdout); 127 args->quit = 1; 128 goto exit; 129 case 'i': args->lines = optarg; break; 130 case 'm': args->output_type = OUTPUT_NODE_MESH; break; 131 case 'n': args->output_type = OUTPUT_NODE_DESCRIPTOR; break; 132 case 'p': args->molparams = optarg; break; 133 case 's': args->lines_in_shtr_format = 1; break; 134 case 'v': args->verbose += (args->verbose < 3); break; 135 case 'w': 136 args->output_type = OUTPUT_NODE_VALUE; 137 res = cstr_to_double(optarg, &args->wavenumber); 138 break; 139 default: res = RES_BAD_ARG; break; 140 } 141 if(res != RES_OK) { 142 if(optarg) { 143 fprintf(stderr, "%s: invalid option argument '%s' -- '%c'\n", 144 argv[0], optarg, opt); 145 } 146 goto error; 147 } 148 } 149 150 #define MANDATORY(Cond, Name, Opt) { \ 151 if(!(Cond)) { \ 152 fprintf(stderr, "%s: %s missing -- option '-%c'\n", argv[0], (Name), (Opt)); \ 153 res = RES_BAD_ARG; \ 154 goto error; \ 155 } \ 156 } (void)0 157 MANDATORY(args->molparams, "molparams", 'p'); 158 MANDATORY(args->lines, "line list", 'i'); 159 #undef MANDATORY 160 161 if(optind < argc) args->tree = argv[optind]; 162 163 exit: 164 return res; 165 error: 166 usage(stderr); 167 goto exit; 168 } 169 170 static void 171 cmd_release(struct cmd* cmd) 172 { 173 ASSERT(cmd); 174 if(cmd->sln) SLN(device_ref_put(cmd->sln)); 175 if(cmd->tree) SLN(tree_ref_put(cmd->tree)); 176 if(cmd->shtr) SHTR(ref_put(cmd->shtr)); 177 if(cmd->molparams) SHTR(isotope_metadata_ref_put(cmd->molparams)); 178 if(cmd->lines) SHTR(line_list_ref_put(cmd->lines)); 179 } 180 181 static res_T 182 load_lines(struct cmd* cmd, const struct args* args) 183 { 184 res_T res = RES_OK; 185 ASSERT(cmd && args); 186 187 if(args->lines_in_shtr_format) { 188 struct shtr_line_list_read_args read_args = SHTR_LINE_LIST_READ_ARGS_NULL; 189 190 /* Loads lines from data serialized by the Star-HITRAN library */ 191 read_args.filename = args->lines; 192 res = shtr_line_list_read(cmd->shtr, &read_args, &cmd->lines); 193 if(res != RES_OK) goto error; 194 195 } else { 196 struct shtr_line_list_load_args load_args = SHTR_LINE_LIST_LOAD_ARGS_NULL; 197 198 /* Loads lines from a file in HITRAN format */ 199 load_args.filename = args->lines; 200 res = shtr_line_list_load(cmd->shtr, &load_args, &cmd->lines); 201 if(res != RES_OK) goto error; 202 } 203 204 exit: 205 return res; 206 error: 207 if(cmd->lines) { SHTR(line_list_ref_put(cmd->lines)); cmd->lines = NULL; } 208 goto exit; 209 } 210 211 static res_T 212 cmd_init(struct cmd* cmd, const struct args* args) 213 { 214 struct sln_device_create_args sln_args = SLN_DEVICE_CREATE_ARGS_DEFAULT; 215 struct sln_tree_read_args tree_args = SLN_TREE_READ_ARGS_NULL; 216 struct shtr_create_args shtr_args = SHTR_CREATE_ARGS_DEFAULT; 217 res_T res = RES_OK; 218 219 ASSERT(cmd && args); 220 221 *cmd = CMD_NULL; 222 223 shtr_args.verbose = args->verbose; 224 res = shtr_create(&shtr_args, &cmd->shtr); 225 if(res != RES_OK) goto error; 226 227 res = shtr_isotope_metadata_load(cmd->shtr, args->molparams, &cmd->molparams); 228 if(res != RES_OK) goto error; 229 230 res = load_lines(cmd, args); 231 if(res != RES_OK) goto error; 232 233 sln_args.verbose = args->verbose; 234 res = sln_device_create(&sln_args, &cmd->sln); 235 if(res != RES_OK) goto error; 236 237 tree_args.metadata = cmd->molparams; 238 tree_args.lines = cmd->lines; 239 if(args->tree) { 240 tree_args.file = NULL; 241 tree_args.filename = args->tree; 242 } else { 243 tree_args.file = stdin; 244 tree_args.filename = "stdin"; 245 } 246 res = sln_tree_read(cmd->sln, &tree_args, &cmd->tree); 247 if(res != RES_OK) goto error; 248 249 cmd->args = *args; 250 251 exit: 252 return res; 253 error: 254 cmd_release(cmd); 255 *cmd = CMD_NULL; 256 goto exit; 257 } 258 259 static res_T 260 print_level_descriptor(const struct cmd* cmd) 261 { 262 /* Stack for visiting the tree depth-first */ 263 struct { 264 const struct sln_node* node; 265 unsigned level; 266 } stack[SLN_TREE_DEPTH_MAX*(SLN_TREE_ARITY_MAX-1/*1st node's child*/)]; 267 int istack = 0; 268 269 /* Node data */ 270 struct sln_node_desc desc = SLN_NODE_DESC_NULL; 271 const struct sln_node* node = NULL; 272 273 /* Level descriptor */ 274 size_t nvertices = 0; 275 size_t nnodes = 0; 276 277 /* Miscellaneous */ 278 unsigned level = 0; 279 res_T res = RES_OK; 280 281 ASSERT(cmd); /* Precondition */ 282 283 /* Push a dummy node which, once pop up, whill mark the end of recursion */ 284 stack[istack].node = NULL; 285 stack[istack].level = UINT_MAX; 286 ++istack; 287 288 node = sln_tree_get_root(cmd->tree); 289 290 while(node) { 291 ASSERT(level <= cmd->args.level); 292 293 if(!sln_node_is_leaf(node) && level < cmd->args.level) { 294 const unsigned nchildren = sln_node_get_child_count(cmd->tree, node); 295 unsigned ichild = 0; 296 297 /* Continue down the tree */ 298 ++level; 299 300 /* Push the node children excepted the 1st */ 301 FOR_EACH(ichild, 1, nchildren) { 302 stack[istack ].node = sln_node_get_child(cmd->tree, node, ichild); 303 stack[istack++].level = level; 304 } 305 306 node = sln_node_get_child(cmd->tree, node, 0); /* Visit the left child */ 307 308 } else { 309 /* The queried level or a leaf is reached, update the descriptor */ 310 if((res = sln_node_get_desc(cmd->tree, node, &desc)) != RES_OK) goto error; 311 nvertices += desc.nvertices; 312 ++nnodes; 313 314 /* Pop the next node */; 315 node = stack[--istack].node; 316 level = stack[ istack].level; 317 } 318 } 319 320 /* Print the level description */ 321 printf("#nodes: %lu\n", (unsigned long)nnodes); 322 printf("#vertices: %lu\n", (unsigned long)nvertices); 323 324 exit: 325 return res; 326 error: 327 goto exit; 328 } 329 330 static const struct sln_node* /* NULL <=> tree is empty */ 331 get_node(const struct cmd* cmd, unsigned* node_depth/*can be NULL*/) 332 { 333 const struct sln_node* node = NULL; 334 unsigned depth = 0; 335 unsigned i = 0; 336 ASSERT(cmd); 337 338 node = sln_tree_get_root(cmd->tree); 339 if(node == NULL) goto exit; /* Tree is empty */ 340 341 FOR_EACH(i, 0, cmd->args.depth) { 342 unsigned nchildren = 0; 343 unsigned ichild = 0; 344 345 if(sln_node_is_leaf(node)) break; 346 347 nchildren = sln_node_get_child_count(cmd->tree, node); 348 ichild = MMIN(cmd->args.descent_path[i], nchildren-1); 349 350 node = sln_node_get_child(cmd->tree, node, ichild); 351 352 ++depth; 353 } 354 355 exit: 356 if(node_depth) *node_depth = depth; 357 return node; 358 } 359 360 static res_T 361 print_node_descriptor(const struct cmd* cmd) 362 { 363 const struct sln_node* node = NULL; 364 struct sln_node_desc desc = SLN_NODE_DESC_NULL; 365 unsigned depth = 0; 366 res_T res = RES_OK; 367 ASSERT(cmd); 368 369 if((node = get_node(cmd, &depth)) == NULL) goto exit; /* tree is empty */ 370 371 res = sln_node_get_desc(cmd->tree, node, &desc); 372 if(res != RES_OK) goto error; 373 374 printf("level: %u\n", depth); 375 printf("#lines: %lu\n", (unsigned long)desc.nlines); 376 printf("#vertices: %lu\n", (unsigned long)desc.nvertices); 377 printf("#children: %u\n", desc.nchildren); 378 379 exit: 380 return res; 381 error: 382 goto exit; 383 } 384 385 static res_T 386 print_mesh(const struct cmd* cmd) 387 { 388 struct sln_mesh mesh = SLN_MESH_NULL; 389 const struct sln_node* node = NULL; 390 size_t i = 0; 391 res_T res = RES_OK; 392 ASSERT(cmd); 393 394 if((node = get_node(cmd, NULL)) == NULL) goto exit; /* tree is empty */ 395 396 res = sln_node_get_mesh(cmd->tree, node, &mesh); 397 if(res != RES_OK) goto error; 398 399 FOR_EACH(i, 0, mesh.nvertices) { 400 printf("%g %g\n", 401 mesh.vertices[i].wavenumber, 402 mesh.vertices[i].ka); 403 } 404 405 exit: 406 return res; 407 error: 408 goto exit; 409 } 410 411 static res_T 412 print_node_value(const struct cmd* cmd) 413 { 414 struct sln_tree_desc tree_desc = SLN_TREE_DESC_NULL; 415 struct sln_mesh mesh = SLN_MESH_NULL; 416 const struct sln_node* node = NULL; 417 double val_mesh = 0; 418 double val_node = 0; 419 res_T res = RES_OK; 420 ASSERT(cmd); 421 422 if((node = get_node(cmd, NULL)) == NULL) goto exit; /* tree is empty */ 423 424 res = sln_node_get_mesh(cmd->tree, node, &mesh); 425 if(res != RES_OK) goto error; 426 427 val_mesh = sln_mesh_eval(&mesh, cmd->args.wavenumber); 428 val_node = sln_node_eval(cmd->tree, node, cmd->args.wavenumber); 429 430 printf("ka(%e) = %e ~ %e\n", cmd->args.wavenumber, val_node, val_mesh); 431 432 res = sln_tree_get_desc(cmd->tree, &tree_desc); 433 if(res != RES_OK) goto error; 434 435 if(tree_desc.mesh_type == SLN_MESH_UPPER && !sln_node_is_leaf(node)) { 436 /* Check that the value of the node is greater than or equal to the sum of 437 * the values of its children */ 438 struct sln_mesh mesh0 = SLN_MESH_NULL; 439 struct sln_mesh mesh1 = SLN_MESH_NULL; 440 const struct sln_node* child0 = NULL; 441 const struct sln_node* child1 = NULL; 442 double val_mesh0 = 0; 443 double val_mesh1 = 0; 444 445 child0 = sln_node_get_child(cmd->tree, node, 0); 446 child1 = sln_node_get_child(cmd->tree, node, 1); 447 if((res = sln_node_get_mesh(cmd->tree, child0, &mesh0)) != RES_OK) goto error; 448 if((res = sln_node_get_mesh(cmd->tree, child1, &mesh1)) != RES_OK) goto error; 449 450 val_mesh0 = sln_mesh_eval(&mesh0, cmd->args.wavenumber); 451 val_mesh1= sln_mesh_eval(&mesh1, cmd->args.wavenumber); 452 453 if(val_mesh < val_mesh0 + val_mesh1) { 454 fprintf(stderr, "error: ka < ka0 + ka1 (ka0=%e; ka1=%e)\n", 455 val_mesh0, val_mesh1); 456 res = RES_BAD_OP; 457 goto error; 458 } 459 } 460 461 exit: 462 return res; 463 error: 464 goto exit; 465 } 466 467 static res_T 468 print_tree_descriptor(const struct cmd* cmd) 469 { 470 struct sln_tree_desc desc = SLN_TREE_DESC_NULL; 471 res_T res = RES_OK; 472 ASSERT(cmd); 473 474 res = sln_tree_get_desc(cmd->tree, &desc); 475 if(res != RES_OK) goto error; 476 477 printf("#lines: %lu\n", (unsigned long)desc.nlines); 478 printf("#nodes: %lu\n", (unsigned long)desc.nnodes); 479 printf("tree depth: %u\n", desc.depth); 480 printf("#vertices: %lu\n", (unsigned long)desc.nvertices); 481 printf("type: %s\n", sln_mesh_type_cstr(desc.mesh_type)); 482 printf("decimation error: %.4e\n", desc.mesh_decimation_err); 483 printf("line profile: %s\n", sln_line_profile_cstr(desc.line_profile)); 484 printf("#lines per leaf: %lu\n", (unsigned long)desc.max_nlines_per_leaf); 485 printf("arity: %u\n", desc.arity); 486 487 exit: 488 return res; 489 error: 490 goto exit; 491 } 492 493 static res_T 494 cmd_run(const struct cmd* cmd) 495 { 496 res_T res = RES_OK; 497 498 switch(cmd->args.output_type) { 499 case OUTPUT_LEVEL_DESCRIPTOR: 500 res = print_level_descriptor(cmd); 501 break; 502 case OUTPUT_NODE_DESCRIPTOR: 503 res = print_node_descriptor(cmd); 504 break; 505 case OUTPUT_NODE_MESH: 506 res = print_mesh(cmd); 507 break; 508 case OUTPUT_NODE_VALUE: 509 res = print_node_value(cmd); 510 break; 511 case OUTPUT_TREE_DESCRIPTOR: 512 res = print_tree_descriptor(cmd); 513 break; 514 default: FATAL("Unreachable code\n"); break; 515 } 516 if(res != RES_OK) goto error; 517 518 exit: 519 return res; 520 error: 521 goto exit; 522 } 523 524 /******************************************************************************* 525 * The program 526 ******************************************************************************/ 527 int 528 main(int argc, char** argv) 529 { 530 struct args args = ARGS_DEFAULT; 531 struct cmd cmd = CMD_NULL; 532 int err = 0; 533 res_T res = RES_OK; 534 535 if((res = args_init(&args, argc, argv)) != RES_OK) goto error; 536 if(args.quit) goto exit; 537 538 if((res = cmd_init(&cmd, &args)) != RES_OK) goto error; 539 if((res = cmd_run(&cmd)) != RES_OK) goto error; 540 541 exit: 542 cmd_release(&cmd); 543 CHK(mem_allocated_size() == 0); 544 return err; 545 error: 546 err = 1; 547 goto exit; 548 }