commit 4219e2e8b5e6afb11d9d2d7ca91b283218b89b64
parent 5d2d3a7d65a5b0c697c9fa622a52fb2e74e13538
Author: Vincent Forest <vincent.forest@meso-star.com>
Date: Wed, 28 Sep 2022 10:42:43 +0200
Start implementating the rngrd_cell_create_phase_fn function
Discrete phase functions are not yet supported
Diffstat:
4 files changed, 243 insertions(+), 2 deletions(-)
diff --git a/README.md b/README.md
@@ -23,6 +23,7 @@ on the [Rad-Net Scattering Function](https://gitlab.com/meso-star/rnsf),
[Star-Buffer](https://gitlab.com/meso-star/star-buffer),
[Star-CorrelatedK](https://gitlab.com/meso-star/star-ck),
[Star-Mesh](https://gitlab.com/meso-star/star-mesh),
+[Star-ScatteringFunctions](https://gitlab.com/meso-star/star-sf),
[Star-UnstructuredVolumetricMesh](https://gitlab.com/meso-star/star-uvm),
[Star-VoXel](https://gitlab.com/meso-star/star-vx)
libraries, and on [OpenMP](https://www.openmp.org) 1.2 to parallelize its
diff --git a/cmake/CMakeLists.txt b/cmake/CMakeLists.txt
@@ -37,6 +37,7 @@ find_package(StarAerosol REQUIRED)
find_package(StarBuffer REQUIRED)
find_package(StarCK REQUIRED)
find_package(StarMesh REQUIRED)
+find_package(StarSF 0.7 REQUIRED)
find_package(StarUVM 0.0 REQUIRED)
find_package(StarVX 0.2 REQUIRED)
@@ -52,6 +53,7 @@ include_directories(
${StarBuffer_INCLUDE_DIR}
${StarCK_INCLUDE_DIR}
${StarMesh_INCLUDE_DIR}
+ ${StarSF_INCLUDE_DIR}
${StarUVM_INCLUDE_DIR}
${StarVX_INCLUDE_DIR})
@@ -90,7 +92,7 @@ rcmake_prepend_path(RNATM_FILES_DOC ${PROJECT_SOURCE_DIR}/../)
add_library(rnatm SHARED ${RNATM_FILES_SRC} ${RNATM_FILES_INC} ${RNATM_FILES_INC_API})
target_link_libraries(rnatm RNSL RNSF RSys StarAerosol StarBuffer StarCK StarMesh
- StarUVM StarVX m)
+ StarSF StarUVM StarVX m)
set_target_properties(rnatm PROPERTIES
COMPILE_FLAGS "${OpenMP_C_FLAGS}"
diff --git a/src/rnatm.h b/src/rnatm.h
@@ -51,6 +51,7 @@
/* Forward declaration of external data types */
struct logger;
struct mem_allocator;
+struct ssf_phase;
enum rnatm_radcoef {
RNATM_RADCOEF_Ka, /* Absorption coefficient */
@@ -115,6 +116,13 @@ struct rnatm_cell_get_radcoef_args {
static const struct rnatm_cell_get_radcoef_args
RNATM_CELL_GET_RADCOEF_ARGS_NULL = RNATM_CELL_GET_RADCOEF_ARGS_NULL__;
+struct rnatm_cell_create_phase_fn_args {
+ struct rnatm_cell cell; /* Cell to query */
+ double barycentric_coords[4]; /* Position into and relative to the cell */
+ double wavelength; /* In nm */
+ double r[2]; /* Random numbers uniformaly distributed in [0, 1[ */
+};
+
struct rnatm_trace_ray_args {
double ray_org[3]; /* Origin of the ray */
double ray_dir[3]; /* Direction of the ray */
@@ -246,6 +254,12 @@ rnatm_cell_get_radcoef
double* k);
RNATM_API res_T
+rnatm_cell_create_phase_fn
+ (struct rnatm* atm,
+ const struct rnatm_cell_create_phase_fn_args* args,
+ struct ssf_phase** phase_fn);
+
+RNATM_API res_T
rnatm_trace_ray
(struct rnatm* rnatm,
const struct rnatm_trace_ray_args* args,
diff --git a/src/rnatm_properties.c b/src/rnatm_properties.c
@@ -28,6 +28,7 @@
#include <star/sars.h>
#include <star/sbuf.h>
#include <star/sck.h>
+#include <star/ssf.h>
#include <rsys/clock_time.h>
#include <rsys/cstr.h>
@@ -42,7 +43,7 @@ check_rnatm_cell_get_radcoef_args
{
double sum_bcoords;
size_t i, n;
- ASSERT(darray_band_size_get(&atm->bands));
+ ASSERT(atm && darray_band_size_get(&atm->bands));
if(!args) return RES_BAD_ARG;
@@ -84,6 +85,49 @@ check_rnatm_cell_get_radcoef_args
return RES_OK;
}
+static res_T
+check_rnatm_cell_create_phase_fn_args
+ (struct rnatm* atm,
+ const struct rnatm_cell_create_phase_fn_args* args)
+{
+ double sum_bcoords;
+ ASSERT(atm);
+
+ if(!args) return RES_BAD_ARG;
+
+ /* Invalid geometric primitive */
+ if(!SUVM_PRIMITIVE_NONE(&args->cell.prim)
+ || args->cell.prim.nvertices != 4) /* Expect a tetraheron */
+ return RES_BAD_ARG;
+
+ /* Invalid component */
+
+ /* Invalid component */
+ if(args->cell.component != RNATM_GAS
+ && args->cell.component >= darray_aerosol_size_get(&atm->aerosols))
+ return RES_BAD_ARG;
+
+ /* Invalid barycentric coordinates */
+ sum_bcoords =
+ args->barycentric_coords[0]
+ + args->barycentric_coords[1]
+ + args->barycentric_coords[2]
+ + args->barycentric_coords[3];
+ if(!eq_eps(sum_bcoords, 1, 1.e-6))
+ return RES_BAD_ARG;
+
+ /* Invalid wavelength */
+ if(args->wavelength < 0)
+ return RES_BAD_ARG;
+
+ /* Invalid random numbers */
+ if(args->r[0] < 0 || args->r[0] >= 1
+ || args->r[1] < 0 || args->r[1] >= 1)
+ return RES_BAD_ARG;
+
+ return RES_OK;
+}
+
static INLINE void
reset_gas_properties(struct gas* gas)
{
@@ -546,6 +590,153 @@ error:
goto exit;
}
+static INLINE res_T
+create_phase_fn_rayleigh(struct rnatm* atm, struct ssf_phase** phase_fn)
+{
+ struct mem_allocator* allocator = NULL;
+ res_T res = RES_OK;
+ ASSERT(atm && phase_fn);
+
+ /* TODO speed up allocation by using a [per thread] fifo allocator */
+ allocator = atm->allocator;
+
+ res = ssf_phase_create(allocator, &ssf_phase_rayleigh, phase_fn);
+ if(res != RES_OK) goto error;
+
+exit:
+ return res;
+error:
+ if(*phase_fn) { SSF(phase_ref_put(*phase_fn)); *phase_fn = NULL; }
+ goto exit;
+}
+
+static INLINE res_T
+create_phase_fn_hg
+ (struct rnatm* atm,
+ const struct rnsf_phase_fn_hg* hg,
+ struct ssf_phase** phase_fn)
+{
+ struct mem_allocator* allocator = NULL;
+ res_T res = RES_OK;
+ ASSERT(atm && hg && phase_fn);
+
+ /* TODO speed up allocation by using a [per thread] fifo allocator */
+ allocator = atm->allocator;
+
+ res = ssf_phase_create(allocator, &ssf_phase_hg, phase_fn);
+ if(res != RES_OK) goto error;
+ res = ssf_phase_hg_setup(*phase_fn, hg->g);
+ if(res != RES_OK) goto error;
+
+exit:
+ return res;
+error:
+ if(*phase_fn) { SSF(phase_ref_put(*phase_fn)); *phase_fn = NULL; }
+ goto exit;
+}
+
+static INLINE res_T
+create_phase_fn_discrete
+ (struct rnatm* atm,
+ const struct rnsf_phase_fn_discrete* discrete,
+ struct ssf_phase** phase_fn)
+{
+ struct mem_allocator* allocator = NULL;
+ res_T res = RES_OK;
+ ASSERT(atm && discrete && phase_fn);
+
+ /* TODO speed up allocation by using a [per thread] fifo allocator */
+ allocator = atm->allocator;
+
+ VFATAL("%s: function not implemented\n", ARG1(FUNC_NAME));
+ (void)allocator;
+ (void)atm;
+ (void)discrete;
+ (void)phase_fn;
+ goto error;
+
+exit:
+ return res;
+error:
+ if(*phase_fn) { SSF(phase_ref_put(*phase_fn)); *phase_fn = NULL; }
+ goto exit;
+}
+
+static res_T
+cell_create_phase_fn_aerosol
+ (struct rnatm* atm,
+ const struct rnatm_cell_create_phase_fn_args* args,
+ struct ssf_phase** phase_fn)
+{
+ struct sbuf_desc sbuf_desc;
+ struct rnsf_phase_fn_hg hg;
+ struct rnsf_phase_fn_discrete discrete;
+ const double* bcoords = NULL;
+ const struct aerosol* aerosol = NULL;
+ const struct rnsf* rnsf = NULL;
+ const struct rnsf_phase_fn* rnsf_phase_fn = NULL;
+ uint32_t phase_fn_id = 0;
+ size_t inode = 0;
+ size_t iphase_fn = 0;
+ int icell_node = 0;
+ res_T res = RES_OK;
+ ASSERT(atm && args && phase_fn);
+ ASSERT(args->cell.component < darray_aerosol_size_get(&atm->aerosols));
+
+ /* Get the aerosol the cell belongs to */
+ aerosol = darray_aerosol_cdata_get(&atm->aerosols) + args->cell.component;
+
+ /* Sample the cell node to consider */
+ bcoords = args->barycentric_coords;
+ if(args->r[0] < bcoords[0]) {
+ icell_node = 0;
+ } else if(args->r[0] < bcoords[0] + bcoords[1]) {
+ icell_node = 1;
+ } else if(args->r[0] < bcoords[0] + bcoords[1] + bcoords[2]) {
+ icell_node = 2;
+ } else {
+ icell_node = 3;
+ }
+
+ /* Retrieve the the phase function id of the node */
+ SBUF(get_desc(aerosol->phase_fn_ids, &sbuf_desc));
+ inode = args->cell.prim.indices[icell_node];
+ ASSERT(inode < sbuf_desc.size);
+ phase_fn_id = *((const uint32_t*)sbuf_desc_at(&sbuf_desc, inode));
+
+ /* Retrieve the spectrally varying phase function of the node */
+ ASSERT(phase_fn_id < darray_phase_fn_size_get(&aerosol->phase_fn_lst));
+ rnsf = darray_phase_fn_cdata_get(&aerosol->phase_fn_lst)[phase_fn_id];
+
+ /* Get the phase function based on the input wavelength */
+ res = rnsf_fetch_phase_fn(rnsf, args->wavelength, args->r[1], &iphase_fn);
+ if(res != RES_OK) goto error;
+ rnsf_phase_fn = rnsf_get_phase_fn(rnsf, iphase_fn);
+
+ /* Create the scattering function */
+ switch(rnsf_phase_fn_get_type(rnsf_phase_fn)) {
+ case RNSF_PHASE_FN_HG:
+ RNSF(phase_fn_get_hg(rnsf_phase_fn, &hg));
+ res = create_phase_fn_hg(atm, &hg, phase_fn);
+ if(res != RES_OK) goto error;
+ break;
+ case RNSF_PHASE_FN_DISCRETE:
+ RNSF(phase_fn_get_discrete(rnsf_phase_fn, &discrete));
+ res = create_phase_fn_discrete(atm, &discrete, phase_fn);
+ if(res != RES_OK) goto error;
+ break;
+ default: FATAL("Unreachable code\n"); break;
+ }
+
+exit:
+ return res;
+error:
+ log_err(atm, "error creating the %s phase function -- %s\n",
+ str_cget(&aerosol->name), res_to_cstr(res));
+ if(*phase_fn) { SSF(phase_ref_put(*phase_fn)); *phase_fn = NULL; }
+ goto exit;
+}
+
/*******************************************************************************
* Exported functions
******************************************************************************/
@@ -582,6 +773,39 @@ error:
goto exit;
}
+res_T
+rnatm_cell_create_phase_fn
+ (struct rnatm* atm,
+ const struct rnatm_cell_create_phase_fn_args* args,
+ struct ssf_phase** out_phase_fn)
+{
+ struct ssf_phase* phase_fn = NULL;
+ res_T res = RES_OK;
+
+ if(!atm || !phase_fn) { res = RES_BAD_ARG; goto error; }
+ res = check_rnatm_cell_create_phase_fn_args(atm, args);
+ if(res != RES_OK) goto error;
+
+ if(args->cell.component != RNATM_GAS) {
+ res = cell_create_phase_fn_aerosol(atm, args, &phase_fn);
+ if(res != RES_OK) goto error;
+ } else {
+ res = create_phase_fn_rayleigh(atm, &phase_fn);
+ if(res != RES_OK) {
+ log_err(atm, "error creating the gas phase function -- %s\n",
+ res_to_cstr(res));
+ goto error;
+ }
+ }
+
+exit:
+ if(out_phase_fn) *out_phase_fn = phase_fn;
+ return res;
+error:
+ if(phase_fn) { SSF(phase_ref_put(phase_fn)); phase_fn = NULL; }
+ goto exit;
+}
+
/*******************************************************************************
* Local functions
******************************************************************************/
