star-hitran

Load line-by-line data from the HITRAN database
git clone git://git.meso-star.fr/star-hitran.git
Log | Files | Refs | README | LICENSE

commit 05e15cd4754a82d08341a4e82a88945b465bc542
parent 8ab70939583d83b146d237aa16e2edcc2d76bb94
Author: Vincent Forest <vincent.forest@meso-star.com>
Date:   Wed,  9 Mar 2022 09:29:16 +0100

Implement shtr_lines_view_<create|ref_put|ref_get>

Diffstat:
Mcmake/CMakeLists.txt | 6++++--
Msrc/shtr.h | 8++++----
Msrc/shtr_lines_list.c | 2+-
Asrc/shtr_lines_view.c | 253+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
4 files changed, 262 insertions(+), 7 deletions(-)

diff --git a/cmake/CMakeLists.txt b/cmake/CMakeLists.txt @@ -48,9 +48,11 @@ set(SHTR_FILES_SRC shtr_log.c shtr_isotope_metadata.c shtr_param.c - shtr_lines_list.c) + shtr_lines_list.c + shtr_lines_view.c) set(SHTR_FILES_INC shtr_c.h + shtr_lines_list_c.h shtr_log.h shtr_param.h) set(SHTR_FILES_INC_API @@ -66,7 +68,7 @@ rcmake_prepend_path(SHTR_FILES_DOC ${PROJECT_SOURCE_DIR}/../) add_library(shtr SHARED ${SHTR_FILES_SRC} ${SHTR_FILES_INC} ${SHTR_FILES_INC_API}) -target_link_libraries(shtr RSys) +target_link_libraries(shtr RSys m) set_target_properties(shtr PROPERTIES DEFINE_SYMBOL SHTR_SHARED_BUILD diff --git a/src/shtr.h b/src/shtr.h @@ -51,10 +51,10 @@ static const struct shtr_create_args SHTR_CREATE_ARGS_DEFAULT = SHTR_CREATE_ARGS_DEFAULT__; struct shtr_lines_view_create_args { - double wavenumbers_range[2]; /* Spectral range */ + double wavenumber_range[2]; /* Spectral range */ /* List of molecule identifiers to consider sorted in ascending order */ - int8_t molecules[SHTR_MAX_MOLECULES_COUNT]; + int32_t molecules[SHTR_MAX_MOLECULES_COUNT]; size_t nmolecules; }; #define SHTR_LINES_VIEW_CREATE_ARGS_NULL__ {{0,0},{0}, 0} @@ -223,9 +223,9 @@ shtr_lines_list_get ******************************************************************************/ SHTR_API res_T shtr_lines_view_create - (const struct shtr_lines_list* lnlst, + (struct shtr_lines_list* lnlst, const struct shtr_lines_view_create_args* args, - struct shtr_lines_view* view); + struct shtr_lines_view** view); SHTR_API res_T shtr_lines_view_ref_get diff --git a/src/shtr_lines_list.c b/src/shtr_lines_list.c @@ -39,7 +39,7 @@ create_lines_list lnlst = MEM_CALLOC(shtr->allocator, 1, sizeof(*lnlst)); if(!lnlst) { - log_err(shtr, "Could not allocate the lnlst data structure.\n"); + log_err(shtr, "Could not allocate the list of lines.\n"); res = RES_MEM_ERR; goto error; } diff --git a/src/shtr_lines_view.c b/src/shtr_lines_view.c @@ -0,0 +1,253 @@ +/* Copyright (C) 2022 CNRS - LMD + * Copyright (C) 2022 |Meso|Star> (contact@meso-star.com) + * Copyright (C) 2022 Université Paul Sabatier - IRIT + * Copyright (C) 2022 Université Paul Sabatier - Laplace + * + * This program is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see <http://www.gnu.org/licenses/>. */ + +#define _POSIX_C_SOURCE 200112L /* nextafter */ + +#include "shtr.h" +#include "shtr_c.h" +#include "shtr_lines_list_c.h" +#include "shtr_log.h" + +#include <rsys/algorithm.h> +#include <rsys/cstr.h> +#include <rsys/dynamic_array_size_t.h> + +#include <math.h> + +struct shtr_lines_view { + struct shtr_lines_list* list; + struct darray_size_t line_ids; /* Indices of the selected lines */ + ref_T ref; +}; + +/******************************************************************************* + * Helper functions + ******************************************************************************/ +static res_T +check_shtr_lines_view_create_args + (struct shtr* shtr, + const struct shtr_lines_view_create_args* args) +{ + if(!args) return RES_BAD_ARG; + if(args->wavenumber_range[0] > args->wavenumber_range[1]) { + log_err(shtr, "Invalid lines view spectral range [%g, %g].\n", + args->wavenumber_range[0], args->wavenumber_range[1]); + return RES_BAD_ARG; + } + return RES_OK; +} + +static int +cmp_line_wavenumber(const void* key, const void* entry) +{ + const double wavenumber = *((double*)key); + const struct shtr_line* line = entry; + if(wavenumber < line->wavenumber) { + return -1; + } else if(wavenumber > line->wavenumber) { + return 1; + } else { + return 0; + } +} + +static res_T +find_line_id_range + (const struct shtr_lines_view* view, + const double wavenumber_range[2], + size_t line_id_range[2]) +{ + const struct shtr_line* lines = NULL; + const struct shtr_line* l0 = NULL; + const struct shtr_line* l1 = NULL; + size_t nlines = 0; + double exclusive_bound = 0; + res_T res = RES_OK; + ASSERT(view && wavenumber_range && line_id_range); + + res = shtr_lines_list_get(view->list, &lines); + if(res != RES_OK) goto error; + res = shtr_lines_list_get_size(view->list, &nlines); + if(res != RES_OK) goto error; + + l0 = search_lower_bound(&wavenumber_range[0], lines, nlines, + sizeof(struct shtr_line), cmp_line_wavenumber); + if(l0) { /* No lower bound is found */ + line_id_range[0] = (size_t)(l0 - lines); + } else { + line_id_range[0] = SIZE_MAX; + line_id_range[1] = 0; + goto exit; + } + + /* The search_lower_bound function returns the first entry that is greater + * than _or_ equal to the search value. This is fine for finding an inclusive + * lower bound, but not for finding an inclusive upper bound like below. In + * the latter case, it is instead necessary to find the first entry greater + * than the value sought and return the previous one. To do this, we + * therefore seek the next representable floating point number directly above + * the upper bound wavenumber.*/ + exclusive_bound = nextafter(wavenumber_range[1], INF); + l1 = search_lower_bound(&exclusive_bound, lines, nlines, + sizeof(struct shtr_line), cmp_line_wavenumber); + if(!l1) { /* No upper bound is found. Take all the lines to the end */ + line_id_range[1] = nlines - 1; + } else { + line_id_range[1] = (size_t)(l1 - lines) - 1; + } + ASSERT(line_id_range[0] <= line_id_range[1]); + ASSERT(line_id_range[1] < nlines); + +exit: + return res; +error: + goto exit; +} + +static int +cmp_int32(const void* a, const void* b) +{ + ASSERT(a && b); + return *((int32_t*)a) - *((int32_t*)b); +} + +/* Return 1 if `molecule' is in `molecules' and 0 otherwise */ +static INLINE int +include_molecule + (const int32_t molecules[], + const size_t nmolecules, + const int32_t molecule) +{ + int32_t* i = NULL; + ASSERT(molecules); + + i = search_lower_bound + (&molecule, molecules, nmolecules, sizeof(int32_t), cmp_int32); + return i != NULL && *i == molecule; +} + +static res_T +select_lines + (struct shtr_lines_view* view, + const struct shtr_lines_view_create_args* args) +{ + const struct shtr_line* lines; + int32_t molecules[SHTR_MAX_MOLECULES_COUNT]; + size_t line_id_range[2]; + size_t i; + res_T res = RES_OK; + ASSERT(view && args); + + /* Nothing to do */ + if(args->nmolecules == 0) goto exit; + + /* Sort in ascending order the indices of the molecules to include */ + memcpy(molecules, args->molecules, args->nmolecules*sizeof(molecules[0])); + qsort(molecules, args->nmolecules, sizeof(molecules[0]), cmp_int32); +#ifndef NDEBUG + FOR_EACH(i, 1, args->nmolecules) { + ASSERT(molecules[i] >= molecules[i-1]); + } +#endif + + res = find_line_id_range(view, args->wavenumber_range, line_id_range); + if(res != RES_OK) goto error; + + lines = darray_line_cdata_get(&view->list->lines); + FOR_EACH(i, line_id_range[0], line_id_range[1]+1) { + if(include_molecule(molecules, args->nmolecules, lines[i].molecule_id)) { + res = darray_size_t_push_back(&view->line_ids, &i); + if(res != RES_OK) { + log_err(view->list->shtr, + "Could not register the line into the view -- %s.\n", + res_to_cstr(res)); + goto error; + } + } + } + +exit: + return res; +error: + darray_size_t_clear(&view->line_ids); + goto exit; +} + +static void +release_lines_view(ref_T* ref) +{ + struct shtr_lines_view* view = CONTAINER_OF(ref, struct shtr_lines_view, ref); + struct shtr_lines_list* list = view->list; + ASSERT(ref); + darray_size_t_release(&view->line_ids); + MEM_RM(list->shtr->allocator, view); + SHTR(lines_list_ref_put(list)); +} + +/******************************************************************************* + * Exported functions + ******************************************************************************/ +res_T +shtr_lines_view_create + (struct shtr_lines_list* list, + const struct shtr_lines_view_create_args* args, + struct shtr_lines_view** out_view) +{ + struct shtr_lines_view* view = NULL; + res_T res = RES_OK; + + if(!list || !out_view) { res = RES_BAD_ARG; goto error; } + res = check_shtr_lines_view_create_args(list->shtr, args); + if(res != RES_OK) goto error; + view = MEM_CALLOC(list->shtr->allocator, 1, sizeof(*view)); + if(!view) { + log_err(list->shtr, "Could not allocate the lines view.\n"); + res = RES_MEM_ERR; + goto error; + } + ref_init(&view->ref); + SHTR(lines_list_ref_get(list)); + view->list = list; + darray_size_t_init(list->shtr->allocator, &view->line_ids); + + res = select_lines(view, args); + if(res != RES_OK) goto error; + +exit: + if(out_view) *out_view = view; + return res; +error: + if(view) { SHTR(lines_view_ref_put(view)); view = NULL; } + goto exit; +} + +res_T +shtr_lines_view_ref_get(struct shtr_lines_view* view) +{ + if(!view) return RES_BAD_ARG; + ref_get(&view->ref); + return RES_OK; +} + +res_T +shtr_lines_view_ref_put(struct shtr_lines_view* view) +{ + if(!view) return RES_BAD_ARG; + ref_put(&view->ref, release_lines_view); + return RES_OK; +}