star-ck

Describe the radiative properties of gas mixtures
git clone git://git.meso-star.fr/star-ck.git
Log | Files | Refs | README | LICENSE

commit 410d6ece4490195046207f48a1f39c0ac2fdad02
parent e87ac331fe392f160e0c71a1669a7e436ad339e0
Author: Vincent Forest <vincent.forest@meso-star.com>
Date:   Wed,  7 Sep 2022 11:55:59 +0200

Add and test the sck_find_bands function

Diffstat:
Msrc/sck.c | 74++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Msrc/sck.h | 9+++++++++
Msrc/test_sck_load.c | 86+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++------
3 files changed, 163 insertions(+), 6 deletions(-)

diff --git a/src/sck.c b/src/sck.c @@ -324,6 +324,23 @@ error: goto exit; } +static INLINE int +cmp_band(const void* key, const void* item) +{ + const struct band* band = item; + double wnum; + ASSERT(key && item); + wnum = *(double*)key; + + if(wnum < band->low) { + return -1; + } else if(wnum >= band->upp) { + return +1; + } else { + return 0; + } +} + static INLINE void hash_quad_pt (struct sha256_ctx* ctx, @@ -527,6 +544,63 @@ error: } res_T +sck_find_bands + (const struct sck* sck, + const double range[2], + size_t ibands[2]) +{ + const struct band* bands = NULL; + const struct band* low = NULL; + const struct band* upp = NULL; + size_t nbands = 0; + res_T res = RES_OK; + + if(!sck || !range || !ibands || range[0] > range[1]) { + res = RES_BAD_ARG; + goto error; + } + + bands = darray_band_cdata_get(&sck->bands); + nbands = darray_band_size_get(&sck->bands); + + low = search_lower_bound(range+0, bands, nbands, sizeof(*bands), cmp_band); + if(low) { + ibands[0] = (size_t)(low - bands); + } else { + /* The submitted range does not overlap any band */ + ibands[0] = SIZE_MAX; + ibands[1] = 0; + goto exit; + } + + if(range[0] == range[1]) goto exit; /* No more to do */ + + upp = search_lower_bound(range+1, bands, nbands, sizeof(*bands), cmp_band); + + /* The submitted range overlaps the remaining bands */ + if(!upp) { + ibands[1] = nbands - 1; + + /* The upper band includes range[1] */ + } else if(upp->low <= range[1]) { + ibands[1] = (size_t)(upp - bands); + + /* The upper band is greater than range[1] and therefre must be rejected */ + } else if(upp->low > range[1]) { + if(upp != bands) { + ibands[1] = (size_t)(upp - bands - 1); + } else { + ibands[0] = SIZE_MAX; + ibands[1] = 0; + } + } + +exit: + return res; +error: + goto exit; +} +res_T sck_band_get_quad_pt (const struct sck_band* sck_band, const size_t iquad_pt, diff --git a/src/sck.h b/src/sck.h @@ -120,6 +120,15 @@ sck_get_band const size_t iband, struct sck_band* band); +/* Returns the range of band indices covered by a given spectral range. The + * returned index range is degenerated (i.e. ibands[0] > ibands[1]) if no band + * is found */ +SCK_API res_T +sck_find_bands + (const struct sck* sck, + const double range[2], /* In nm. Limits are inclusive */ + size_t ibands[2]); /* Range of overlaped bands. Limits are inclusive */ + SCK_API res_T sck_band_get_quad_pt (const struct sck_band* band, diff --git a/src/test_sck_load.c b/src/test_sck_load.c @@ -14,6 +14,8 @@ * 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 /* nextarfter */ + #include "sck.h" #include <rsys/math.h> @@ -66,8 +68,8 @@ check_sck_load CHK(sck_band_sample_quad_pt(&band, 0, NULL) == RES_BAD_ARG); FOR_EACH(iband, 0, nbands) { - const double low = (double)iband; - const double upp = (double)(iband+1); + const double low = (double)(iband+1); + const double upp = (double)(iband+2); const size_t nqpts = iband + 1; size_t iqpt; double sum[NMOMENTS] = {0}; @@ -157,8 +159,8 @@ write_sck CHK(fwrite(&nnodes, sizeof(nnodes), 1, fp) == 1); FOR_EACH(iband, 0, nbands) { - const double low = (double)iband; - const double upp = (double)(iband+1); + const double low = (double)(iband+1); + const double upp = (double)(iband+2); const uint64_t nqpts = iband + 1; uint64_t iqpt; @@ -413,6 +415,78 @@ test_load_files(struct sck* sck, int argc, char** argv) } } +static void +test_find(struct sck* sck) +{ + size_t ibands[2]; + double range[2]; + FILE* fp; + + CHK(fp = tmpfile()); + write_sck(fp, 4096, 10, 1); + rewind(fp); + CHK(sck_load_stream(sck, fp, NULL) == RES_OK); + + range[0] = 0; + range[1] = 10; + CHK(sck_find_bands(NULL, range, ibands) == RES_BAD_ARG); + CHK(sck_find_bands(sck, NULL, ibands) == RES_BAD_ARG); + CHK(sck_find_bands(sck, range, NULL) == RES_BAD_ARG); + CHK(sck_find_bands(sck, range, ibands) == RES_OK); + CHK(ibands[0] == 0 && ibands[1] == 9); + + range[0] = 10; + range[1] = 0; + CHK(sck_find_bands(sck, range, ibands) == RES_BAD_ARG); + + range[0] = 11; + range[1] = 12; + CHK(sck_find_bands(sck, range, ibands) == RES_OK); + CHK(ibands[0] > ibands[1]); + + range[0] = 0; + range[1] = nextafter(1, 0); + CHK(sck_find_bands(sck, range, ibands) == RES_OK); + CHK(ibands[0] > ibands[1]); + + range[0] = 1; + range[1] = 1; + CHK(sck_find_bands(sck, range, ibands) == RES_OK); + CHK(ibands[0] == 0 && ibands[1] == 0); + + range[0] = 0; + range[1] = 1; + CHK(sck_find_bands(sck, range, ibands) == RES_OK); + CHK(ibands[0] == 0 && ibands[1] == 0); + + range[0] = 1; + range[1] = nextafterf(2, 1); + CHK(sck_find_bands(sck, range, ibands) == RES_OK); + CHK(ibands[0] == 0 && ibands[1] == 0); + + range[0] = 2; + range[1] = 20; + CHK(sck_find_bands(sck, range, ibands) == RES_OK); + CHK(ibands[0] == 1 && ibands[1] == 9); + + range[0] = 1.5; + range[1] = 2; + CHK(sck_find_bands(sck, range, ibands) == RES_OK); + CHK(ibands[0] == 0 && ibands[1] == 1); + + range[0] = 3.1; + range[1] = nextafter(6, 0); + CHK(sck_find_bands(sck, range, ibands) == RES_OK); + CHK(ibands[0] == 2 && ibands[1] == 4); + + range[0] = 3.1; + range[1] = 7; + CHK(sck_find_bands(sck, range, ibands) == RES_OK); + CHK(ibands[0] == 2 && ibands[1] == 6); + + CHK(fclose(fp) == 0); +} + /******************************************************************************* * Main function ******************************************************************************/ @@ -425,12 +499,12 @@ main(int argc, char** argv) args.verbose = 1; CHK(sck_create(&args, &sck) == RES_OK); - - if(argc > 1) { +if(argc > 1) { test_load_files(sck, argc, argv); } else { test_load(sck); test_load_fail(sck); + test_find(sck); } CHK(sck_ref_put(sck) == RES_OK);