htsky

Load and structure a vertically stratified atmosphere
git clone git://git.meso-star.fr/htsky.git
Log | Files | Refs | README | LICENSE

commit 1dfdf41d8f1e3467bc8dbfc16e46ecc4b033f4a1
parent fbe01d03cfdb21b125de724ce1acda069507e05c
Author: Vincent Forest <vincent.forest@meso-star.com>
Date:   Mon, 13 Nov 2023 12:18:25 +0100

Merge branch 'release_0.3.1'

Diffstat:
MREADME.md | 9+++++++++
Mconfig.mk | 2+-
Msrc/htsky.c | 4++++
Msrc/htsky_cloud.c | 47++++++++++++++++++++++++++++++++++++++---------
4 files changed, 52 insertions(+), 10 deletions(-)

diff --git a/README.md b/README.md @@ -23,6 +23,15 @@ Edit config.mk as needed, then run: ## Release notes +### Version 0.3.1 + +- Fix invalid memory read for clouds with irregular Z dimensions. This + bug was triggered by clouds that were taller than they were wide, i.e. + with the largest Z dimension. It caused a crash when accessing an + out-of-bounds cell in the data structure used to speed up cloud + indexing. In other situations, i.e. when clouds were wider than they + were high, it had no impact on calculations. + ### Version 0.3 - Replace CMake by Makefile as build system. diff --git a/config.mk b/config.mk @@ -1,4 +1,4 @@ -VERSION = 0.3.0 +VERSION = 0.3.1 PREFIX = /usr/local LIB_TYPE = SHARED diff --git a/src/htsky.c b/src/htsky.c @@ -829,8 +829,10 @@ htsky_fetch_raw_property * position in the svx2htcp_z Look Up Table and use the LUT to define the * voxel index into the HTCP descriptor */ const struct split* splits = darray_split_cdata_get(&sky->svx2htcp_z); + const size_t nsplits = darray_split_size_get(&sky->svx2htcp_z); const size_t ilut = (size_t) ((pos_cs[2] - cloud_desc->lower[2]) / sky->lut_cell_sz); + if(ilut >= nsplits) FATAL("LUT index out of range\n"); ivox[2] = splits[ilut].index + (pos_cs[2] > splits[ilut].height); } @@ -935,7 +937,9 @@ htsky_fetch_temperature(const struct htsky* sky, const double pos[3]) * position in the svx2htcp_z Look Up Table and use the LUT to define the * voxel index into the HTCP descripptor */ const struct split* splits = darray_split_cdata_get(&sky->svx2htcp_z); + const size_t nsplits = darray_split_size_get(&sky->svx2htcp_z); const size_t ilut = (size_t)((pos_cs[2] - desc->lower[2]) / sky->lut_cell_sz); + if(ilut >= nsplits) FATAL("LUT index out of range\n"); ivox[2] = splits[ilut].index + (pos_cs[2] > splits[ilut].height); } diff --git a/src/htsky_cloud.c b/src/htsky_cloud.c @@ -214,6 +214,7 @@ compute_k_bounds_irregular_z double Ca_avg; double Cs_avg; double pos_z; + double lut_low, lut_upp; size_t ilut_low, ilut_upp; size_t ilut; size_t ivox_z_prev = SIZE_MAX; @@ -241,8 +242,20 @@ compute_k_bounds_irregular_z /* Compute the *inclusive* bounds of the indices of the LUT cells * overlapped by the SVX voxel */ - ilut_low = (size_t)((vox_svx_low[2]-sky->htcp_desc.lower[2])/sky->lut_cell_sz); - ilut_upp = (size_t)((vox_svx_upp[2]-sky->htcp_desc.lower[2])/sky->lut_cell_sz); + lut_low = (vox_svx_low[2]-sky->htcp_desc.lower[2])/sky->lut_cell_sz; + lut_upp = (vox_svx_upp[2]-sky->htcp_desc.lower[2])/sky->lut_cell_sz; + ilut_low = (size_t)lut_low; + ilut_upp = (size_t)lut_upp; + + /* A LUT coordinate equal to its search index means that this coordinate lies + * strictly on the boundary of a LUT cell. If this coordinate is the upper + * limit of a coordinate range, we subtract one from the LUT index to ensure + * that the corresponding cell is indeed included in the submitted range + * coordinates */ + if(lut_upp == ilut_upp) { + --ilut_upp; + } + ASSERT(ilut_low <= ilut_upp); /* Prepare the iteration over the LES voxels overlapped by the SVX voxel */ @@ -250,8 +263,10 @@ compute_k_bounds_irregular_z ka_bounds[1] = ks_bounds[1] = kext_bounds[1] =-DBL_MAX; ivox_z_prev = SIZE_MAX; pos_z = vox_svx_low[2]; - ASSERT(pos_z >= (double)ilut_low * sky->lut_cell_sz); - ASSERT(pos_z <= (double)ilut_upp * sky->lut_cell_sz); + ASSERT(vox_svx_low[2] >= + sky->htcp_desc.lower[2] + sky->lut_cell_sz*(double)(ilut_low)); + ASSERT(vox_svx_upp[2] <= + sky->htcp_desc.lower[2] + sky->lut_cell_sz*(double)(ilut_upp+1)); /* Iterate over the LUT cells overlapped by the voxel */ FOR_EACH(ilut, ilut_low, ilut_upp+1) { @@ -414,6 +429,7 @@ compute_xh2o_range_irregular_z size_t ilut_low, ilut_upp; size_t ilut; size_t ivox_z_prev; + double lut_low, lut_upp; double low[2], upp[2]; double pos_z; double ipart; @@ -435,17 +451,30 @@ compute_xh2o_range_irregular_z /* Compute the *inclusive* bounds of the indices of the LUT cells * overlapped by the SVX voxel */ - ilut_low = (size_t)((vox_svx_low[2]-sky->htcp_desc.lower[2])/sky->lut_cell_sz); - ilut_upp = (size_t)((vox_svx_upp[2]-sky->htcp_desc.lower[2])/sky->lut_cell_sz); + lut_low = (vox_svx_low[2]-sky->htcp_desc.lower[2])/sky->lut_cell_sz; + lut_upp = (vox_svx_upp[2]-sky->htcp_desc.lower[2])/sky->lut_cell_sz; + ilut_low = (size_t)lut_low; + ilut_upp = (size_t)lut_upp; + + /* A LUT coordinate equal to its search index means that this coordinate lies + * strictly on the boundary of a LUT cell. If this coordinate is the upper + * limit of a coordinate range, we subtract one from the LUT index to ensure + * that the corresponding cell is indeed included in the submitted range + * coordinates */ + if(lut_upp == ilut_upp) { + --ilut_upp; + } + ASSERT(ilut_low <= ilut_upp); /* Prepare the iteration over the LES voxels overlapped by the SVX voxel */ x_h2o_range[0] = DBL_MAX; x_h2o_range[1] =-DBL_MAX; ivox_z_prev = SIZE_MAX; - pos_z = vox_svx_low[2]; - ASSERT(pos_z >= (double)ilut_low * sky->lut_cell_sz); - ASSERT(pos_z <= (double)ilut_upp * sky->lut_cell_sz); + ASSERT(vox_svx_low[2] >= + sky->htcp_desc.lower[2] + sky->lut_cell_sz*(double)(ilut_low)); + ASSERT(vox_svx_upp[2] <= + sky->htcp_desc.lower[2] + sky->lut_cell_sz*(double)(ilut_upp+1)); /* Iterate over the LUT cells overlapped by the voxel */ FOR_EACH(ilut, ilut_low, ilut_upp+1) {