commit d04b25ded57001414af7eb4a3f863000b8460add
parent df29d4ea741ce6653129dfc6b32d9e955cd2e268
Author: Vincent Forest <vincent.forest@meso-star.com>
Date: Fri, 25 Sep 2015 11:13:50 +0200
Integrate the Gebhart Factors for each spectral band
Add the spectral_band identifier to the getter of the estimator status.
Diffstat:
3 files changed, 60 insertions(+), 34 deletions(-)
diff --git a/src/sgf.h b/src/sgf.h
@@ -110,14 +110,17 @@ SGF_API res_T
sgf_device_ref_put
(struct sgf_device* dev);
-/* TODO comment it */
+/* Numerical estimation of the Gebhart Factor between `primitive_id' and the
+ * other scene primitives. An active S3D_TRACE/S3D_GET_PRIMITIVE session must
+ * be active on `desc->scene'. This function is thread safe and can be invoked
+ * by several threads. */
SGF_API res_T
sgf_integrate
(struct sgf_device* dev,
struct ssp_rng* rng,
const size_t primitive_id,
struct sgf_scene_desc* desc,
- const size_t steps_count,
+ const size_t steps_count, /* #realisations */
struct sgf_estimator** out_estimator);
SGF_API res_T
@@ -132,6 +135,7 @@ SGF_API res_T
sgf_estimator_get_status
(struct sgf_estimator* estimator,
const size_t primitive_id,
+ const size_t spectral_band_id,
struct sgf_status* status);
END_DECLS
diff --git a/src/sgf_estimator.c b/src/sgf_estimator.c
@@ -49,7 +49,9 @@ struct accum {
/* Estimator of the Gebhart Factor of a given face to all the other ones */
struct sgf_estimator {
struct darray_accum buf; /* Integrated radiative flux of each primitive */
- size_t nsteps; /* Number of realisations */
+ size_t nsteps; /* # realisations */
+ size_t nprims; /* # primitives */
+ size_t nbands; /* # spectral bands */
struct sgf_device* dev;
ref_T ref;
@@ -66,22 +68,17 @@ check_scene_desc(struct sgf_scene_desc* desc)
}
static FINLINE void
-accum_weight
- (struct sgf_estimator* estimator,
- const size_t iface,
- const double weight)
+accum_weight(struct accum* accums, const size_t iface, const double weight)
{
- struct accum* acc;
- ASSERT(estimator);
- acc = darray_accum_data_get(&estimator->buf) + iface;
- acc->radiative_flux += weight;
- acc->sqr_radiative_flux += weight * weight;
+ ASSERT(accums);
+ accums[iface].radiative_flux += weight;
+ accums[iface].sqr_radiative_flux += weight * weight;
}
/* Return RES_BAD_OP on numerical issue. i.e. the radiative path is skipped */
static res_T
gebhart_radiative_path
- (struct sgf_estimator* estimator,
+ (struct accum* accums,
struct ssp_rng* rng,
const size_t primitive_id,
struct sgf_scene_desc* desc)
@@ -159,7 +156,7 @@ gebhart_radiative_path
if(transmissivity > trans_min) {
const double weight = transmissivity * emissivity;
- accum_weight(estimator, prim.scene_prim_id, weight);
+ accum_weight(accums, prim.scene_prim_id, weight);
transmissivity = transmissivity * (1.0 - emissivity);
#ifndef NDEBUG
sum_radiative_flux += weight;
@@ -168,7 +165,7 @@ gebhart_radiative_path
/* Russian roulette */
if(ssp_rng_canonical(rng) < emissivity) {
const double weight = transmissivity;
- accum_weight(estimator, prim.scene_prim_id, weight);
+ accum_weight(accums, prim.scene_prim_id, weight);
#ifndef NDEBUG
sum_radiative_flux += weight;
#endif
@@ -176,8 +173,7 @@ gebhart_radiative_path
}
}
- if(reflectivity <= 0.0)
- break;
+ if(reflectivity <= 0.0) break;
proba_reflec_spec = specularity / reflectivity;
f3_normalize(normal, hit.normal);
@@ -259,7 +255,8 @@ sgf_integrate
struct sgf_estimator** out_estimator)
{
struct sgf_estimator* estimator = NULL;
- size_t istep = 0;
+ size_t istep;
+ size_t ispectral_band;
size_t nprims;
int mask;
res_T res = RES_OK;
@@ -307,12 +304,18 @@ sgf_integrate
memset(darray_accum_data_get(&estimator->buf), 0,
darray_accum_size_get(&estimator->buf)*sizeof(struct accum));
- /* Invoke the Gebhart factor integration. TODO interate per spectral bands */
- FOR_EACH(istep, 0, steps_count) {
- res = gebhart_radiative_path(estimator, rng, primitive_id, desc);
- if(res != RES_OK) goto error;
+ /* Invoke the Gebhart factor integration. */
+ FOR_EACH(ispectral_band, 0, desc->spectral_bands_count) {
+ struct accum* accums =
+ darray_accum_data_get(&estimator->buf) + ispectral_band * nprims;
+ FOR_EACH(istep, 0, steps_count) {
+ res = gebhart_radiative_path(accums, rng, primitive_id, desc);
+ if(res != RES_OK) goto error;
+ }
}
estimator->nsteps = steps_count;
+ estimator->nprims = nprims;
+ estimator->nbands = desc->spectral_bands_count;
exit:
if(out_estimator) *out_estimator = estimator;
@@ -342,23 +345,34 @@ res_T
sgf_estimator_get_status
(struct sgf_estimator* estimator,
const size_t iprimitive,
+ const size_t ispectral_band,
struct sgf_status* status)
{
struct accum* acc;
size_t nsteps;
+ size_t iacc;
if(!estimator || !status)
return RES_BAD_ARG;
- if(iprimitive >= darray_accum_size_get(&estimator->buf)) {
+ if(iprimitive >= estimator->nprims) {
+ if(estimator->dev->verbose) {
+ logger_print(estimator->dev->logger, LOG_ERROR,
+ "Out of bound primitive index `%lu'\n", iprimitive);
+ }
+ return RES_BAD_ARG;
+ }
+ if(ispectral_band >= estimator->nbands) {
if(estimator->dev->verbose) {
logger_print(estimator->dev->logger, LOG_ERROR,
- "Invalid primitive index `%lu'\n", iprimitive);
+ "Out of bound spectral band index `%lu'\n", ispectral_band);
}
return RES_BAD_ARG;
}
- acc = darray_accum_data_get(&estimator->buf) + iprimitive;
+ iacc = ispectral_band * estimator->nprims + iprimitive;
+ ASSERT(iacc < darray_accum_size_get(&estimator->buf));
+ acc = darray_accum_data_get(&estimator->buf) + iacc;
nsteps = estimator->nsteps;
status->E = acc->radiative_flux / (double)nsteps; /* Expected value */
diff --git a/src/test_sgf_cube.c b/src/test_sgf_cube.c
@@ -216,14 +216,22 @@ main(int argc, char** argv)
CHECK(sgf_integrate(NULL, rngs[0], 0, &scn_desc, NSTEPS, &estimator), RES_BAD_ARG);
CHECK(sgf_integrate(sgf, rngs[0], 0, &scn_desc, NSTEPS, &estimator), RES_OK);
- CHECK(sgf_estimator_get_status(NULL, SIZE_MAX, NULL), RES_BAD_ARG);
- CHECK(sgf_estimator_get_status(estimator, SIZE_MAX, NULL), RES_BAD_ARG);
- CHECK(sgf_estimator_get_status(NULL, 0, NULL), RES_BAD_ARG);
- CHECK(sgf_estimator_get_status(estimator, 0, NULL), RES_BAD_ARG);
- CHECK(sgf_estimator_get_status(NULL, SIZE_MAX, status), RES_BAD_ARG);
- CHECK(sgf_estimator_get_status(estimator, SIZE_MAX, status), RES_BAD_ARG);
- CHECK(sgf_estimator_get_status(NULL, 0, status), RES_BAD_ARG);
- CHECK(sgf_estimator_get_status(estimator, 0, status), RES_OK);
+ CHECK(sgf_estimator_get_status(NULL, SIZE_MAX, SIZE_MAX, NULL), RES_BAD_ARG);
+ CHECK(sgf_estimator_get_status(estimator, SIZE_MAX, SIZE_MAX, NULL), RES_BAD_ARG);
+ CHECK(sgf_estimator_get_status(NULL, 0, SIZE_MAX, NULL), RES_BAD_ARG);
+ CHECK(sgf_estimator_get_status(estimator, 0, SIZE_MAX, NULL), RES_BAD_ARG);
+ CHECK(sgf_estimator_get_status(NULL, SIZE_MAX, SIZE_MAX, status), RES_BAD_ARG);
+ CHECK(sgf_estimator_get_status(estimator, SIZE_MAX, SIZE_MAX, status), RES_BAD_ARG);
+ CHECK(sgf_estimator_get_status(NULL, 0, SIZE_MAX, status), RES_BAD_ARG);
+ CHECK(sgf_estimator_get_status(estimator, 0, SIZE_MAX, status), RES_BAD_ARG);
+ CHECK(sgf_estimator_get_status(NULL, SIZE_MAX, 0, NULL), RES_BAD_ARG);
+ CHECK(sgf_estimator_get_status(estimator, SIZE_MAX, 0, NULL), RES_BAD_ARG);
+ CHECK(sgf_estimator_get_status(NULL, 0, 0, NULL), RES_BAD_ARG);
+ CHECK(sgf_estimator_get_status(estimator, 0, 0, NULL), RES_BAD_ARG);
+ CHECK(sgf_estimator_get_status(NULL, SIZE_MAX, 0, status), RES_BAD_ARG);
+ CHECK(sgf_estimator_get_status(estimator, SIZE_MAX, 0, status), RES_BAD_ARG);
+ CHECK(sgf_estimator_get_status(NULL, 0, 0, status), RES_BAD_ARG);
+ CHECK(sgf_estimator_get_status(estimator, 0, 0, status), RES_OK);
CHECK(status->nsteps, NSTEPS);
@@ -244,7 +252,7 @@ main(int argc, char** argv)
CHECK(sgf_integrate(sgf, rngs[ithread], iprim, &scn_desc, NSTEPS, &est), RES_OK);
FOR_EACH(iprim2, 0, nprims) {
- CHECK(sgf_estimator_get_status(est, iprim2, row + iprim2), RES_OK);
+ CHECK(sgf_estimator_get_status(est, iprim2, 0, row + iprim2), RES_OK);
CHECK(row[iprim2].nsteps, NSTEPS);
}
CHECK(sgf_estimator_ref_put(est), RES_OK);
