commit 682d283e91c2c142205d9eda4e7f33597fad08aa
parent 2a7b25faa69b6261e2a04b78ec8870f731974a87
Author: Vincent Forest <vincent.forest@meso-star.com>
Date: Tue, 12 Apr 2016 15:57:57 +0200
Clean-up the 4V/S realization code
Fix possible memory leaks when an error occurs. Refactor some code.
Diffstat:
| M | src/realization.c | | | 136 | +++++++++++++++++++++++++++++++++++-------------------------------------------- |
| M | src/realization.h | | | 12 | +++++------- |
| M | src/s4vs.c | | | 211 | ++++++++++++++++++++++++++++++++++++++++--------------------------------------- |
3 files changed, 173 insertions(+), 186 deletions(-)
diff --git a/src/realization.c b/src/realization.c
@@ -38,103 +38,87 @@
#define MAX_FAILURES 10
void
-realization(void* length, struct ssp_rng* rng, void* context)
+realization(void* out_length, struct ssp_rng* rng, void* context)
{
struct context* ctx = (struct context*)context;
struct s3d_attrib attrib;
struct s3d_primitive prim;
- double lambda;
- float normal[3], direction[4], origin[3], range[2], st[2];
+ float normal[3], dir[4], pos[3], range[2], st[2];
struct s3d_hit hit;
+ double length = 0;
+ double lambda = 0;
int nfailures = 0;
- SMC_DOUBLE(length) = 0;
-
- /* sample a starting point on the scene surface to start the random walk from */
- for(;;) {
- /* sample the scene */
+ do { /* Sample the scene surface to define the random walk starting point */
float u, v, w;
- do { u = (float)ssp_rng_canonical(rng); } while(u >= 1);
- do { v = (float)ssp_rng_canonical(rng); } while(v >= 1);
- do { w = (float)ssp_rng_canonical(rng); } while(w >= 1);
- /* get a location: primitive ID and parametric coordinates */
+
+ /* Sample a surface location, i.e. primitive ID and parametric coordinates */
+ u = ssp_rng_canonical_float(rng);
+ v = ssp_rng_canonical_float(rng);
+ w = ssp_rng_canonical_float(rng);
S3D(scene_sample(ctx->scene, u, v, w, &prim, st));
/* retrieve the sampled geometric normal and position */
S3D(primitive_get_attrib(&prim, S3D_GEOMETRY_NORMAL, st, &attrib));
f3_normalize(normal, attrib.value);
S3D(primitive_get_attrib(&prim, S3D_POSITION, st, &attrib));
- f3_set(origin, attrib.value);
-
- /* cosine weighted sampling of the hemisphere around the sampled normal */
- ssp_ran_hemisphere_cos(rng, normal, direction);
-
- /* trace a ray from the sampled location */
- range[0] = 1.e-6f; range[1] = FLT_MAX;
- for(;;) {
- /* find first intersection with geometry to determine its distance */
- S3D(scene_trace_ray(ctx->scene, origin, direction, range, &hit));
- /* handle auto-intersection: hitting the same primitive means auto-intersection */
- if(S3D_PRIMITIVE_EQ(&hit.prim, &prim)) {
- /* just increase the minimum distance and retry */
- range[0] += 1.e-6f;
- }
- else break; /* standard case: OK */
- }
+ f3_set(pos, attrib.value);
- f3_normalize(hit.normal, hit.normal);
+ /* Cosine weighted sampling of the hemisphere around the sampled normal */
+ ssp_ran_hemisphere_cos(rng, normal, dir);
- /* handle precision issues (no intersection means geometry leakage: error) */
- if(S3D_HIT_NONE(&hit)) {
- if(++nfailures >= MAX_FAILURES) {
- /* too many errors: stop */
- ctx->exit_failure = 1;
- break;
- }
- }
- else break; /* standard case: OK */
- }
- if(ctx->exit_failure)
- return;
+ /* Find the 1st hit from the sampled location along the sampled direction */
+ range[0] = 1.e-6f;
+ range[1] = FLT_MAX;
+ do {
+ S3D(scene_trace_ray(ctx->scene, pos, dir, range, &hit));
+ range[0] += 1.e-6f;
+ } while(S3D_PRIMITIVE_EQ(&hit.prim, &prim)); /* Handle self-hit */
+
+ /* No intersection <=> numerical imprecision or geometry leakage */
+ } while(S3D_HIT_NONE(&hit) && ++nfailures < MAX_FAILURES);
- /* here the starting point and a propagation direction have been sampled */
+ /* Too many ray-tracing failures => the geometry is inconsistent */
+ if(nfailures >= MAX_FAILURES) goto error;
+
+ /* Here the starting point and a propagation direction have been sampled */
/* and we have determined the distance of the geometry in this direction */
- /* start diffuse random walk */
- for(;;) {
- /* sample a free path according to ks */
- lambda = ssp_ran_exp(rng, ctx->ks);
+ while(1) { /* Here we go for the diffuse random walk */
+ lambda = ssp_ran_exp(rng, ctx->ks); /* Sample a length according to ks */
+
if(lambda >= hit.distance) {
- /* lambda exceeds the geometry distance: the random walk ends on geometry */
- SMC_DOUBLE(length) += hit.distance;
- break; /* the standard way of leaving the loop */
- } else {
- /* lambda doesn't exceed the geometry distance */
- /* the random walk can continue after a diffusion */
- SMC_DOUBLE(length) += lambda;
- /* new location after lambda */
- f3_add(origin, origin, f3_mulf(direction, direction, (float)lambda));
-
- do {
- /* sample a new direction */
- ssp_ran_sphere_uniform(rng, direction);
- range[0] = 0;
- /* find the first intersection with the geometry */
- /* to refresh the geometry's distance */
- S3D(scene_trace_ray(ctx->scene, origin, direction, range, &hit));
-
- /* handle precision issues (no intersection means geometry leakage: error) */
- if(S3D_HIT_NONE(&hit)) {
- if(++nfailures >= MAX_FAILURES) {
- /* too many errors: stop */
- ctx->exit_failure = 1;
- break;
- }
- }
- /* should be an infinite loop, or there is a geometry leakage */
- } while(S3D_HIT_NONE(&hit));
+ /* Lambda exceeds the geometry distance: stop the random walk */
+ length += hit.distance;
+ break;
}
+
+ /* lambda doesn't exceed the geometry distance */
+ /* the random walk can continue after a diffusion */
+ length += lambda;
+ f3_add(pos, pos, f3_mulf(dir, dir, (float)lambda));
+
+ do {
+ /* sample a new direction */
+ ssp_ran_sphere_uniform(rng, dir);
+ range[0] = 0;
+ /* find the first intersection with the geometry */
+ /* to refresh the geometry's distance */
+ S3D(scene_trace_ray(ctx->scene, pos, dir, range, &hit));
+
+ /* No intersection <=> numerical imprecision or geometry leakage */
+ } while(S3D_HIT_NONE(&hit) && ++nfailures < MAX_FAILURES);
+
+ /* Too many ray-tracing failures => the geometry is inconsistent */
+ if(nfailures >= MAX_FAILURES) goto error;
}
- ATOMIC_ADD(&ctx->nfailures, nfailures);
+exit:
+ ATOMIC_ADD(&ctx->nfailures, nfailures); /* Register the number of failures */
+ SMC_DOUBLE(out_length) = length;
+ return;
+error:
+ ctx->error = 1; /* Notify that an error occurs */
+ length = 0; /* Reset the path length */
+ goto exit;
}
diff --git a/src/realization.h b/src/realization.h
@@ -26,8 +26,8 @@
* The fact that you are presently reading this means that you have had
* knowledge of the CeCILL license and that you accept its terms. */
-#ifndef REALIZATION_H_
-#define REALIZATION_H_
+#ifndef REALIZATION_H
+#define REALIZATION_H
/* forward definition */
struct ssp_rng;
@@ -37,14 +37,12 @@ struct context {
struct s3d_scene* scene;
double ks;
int nfailures;
- char exit_failure;
+ int error;
};
/*******************************************************************************
* MC realization function
******************************************************************************/
+extern void realization(void* length, struct ssp_rng* rng, void* context);
-void
-realization(void* length, struct ssp_rng* rng, void* context);
-
-#endif
+#endif /* REALIZATION_H */
diff --git a/src/s4vs.c b/src/s4vs.c
@@ -27,6 +27,7 @@
* knowledge of the CeCILL license and that you accept its terms. */
#include <rsys/clock_time.h>
+#include <rsys/cstr.h>
#include <rsys/float3.h>
#include <rsys/mem_allocator.h>
@@ -42,174 +43,178 @@ compute_4v_s(struct s3d_scene* scene, const size_t max_steps, const double ks)
char buf[512];
struct time begin, end, elapsed;
struct context ctx;
- float S, V, expected;
- struct smc_device* smc;
+ float S, V, reference;
+ struct smc_device* smc = NULL;
struct smc_integrator integrator;
struct smc_estimator* estimator = NULL;
struct smc_estimator_status estimator_status;
double length;
double sig_length;
+ int mask;
res_T res = RES_OK;
-
- ASSERT(max_steps > 0 && ks > 0);
+ ASSERT(scene && max_steps > 0 && ks > 0);
S3D(scene_begin_session(scene, S3D_SAMPLE|S3D_TRACE));
- /* compute the expected result using a mesh-based method */
+ /* Compute the expected result using a mesh-based method */
S3D(scene_compute_area(scene, &S));
S3D(scene_compute_volume(scene, &V));
if(eq_epsf(S, 0, 1.e-6f) || S < 0) {
- printf("No surface to sample. Is the scene empty?\n");
+ fprintf(stderr, "No surface to sample. Is the scene empty?\n");
+ res = RES_BAD_ARG;
goto error;
}
if(eq_epsf(V, 0, 1.e-6f) || V < 0) {
- printf("Invalid volume \"%.2f\". The scene might not match the prerequisites:\n"
+ fprintf(stderr,
+ "Invalid volume \"%.2f\". The scene might not match the prerequisites:\n"
"it must be closed and its normals must point *into* the volume.\n", V);
+ res = RES_BAD_ARG;
goto error;
}
- expected = 4 * V / S;
+ reference = 4*V/S;
- /* initialize context for MC computation */
+ /* Initialize context for MC computation */
ctx.scene = scene;
ctx.nfailures = 0;
- ctx.exit_failure = 0;
+ ctx.error = 0;
ctx.ks = ks;
- /* setup Star-MC */
- integrator.integrand = &realization; /* realization function */
- integrator.type = &smc_double; /* realization type */
- integrator.max_steps = max_steps; /* realization count */
+ /* Setup Star-MC */
+ SMC(device_create(NULL, NULL, SMC_NTHREADS_DEFAULT, NULL, &smc));
+ integrator.integrand = &realization; /* Realization function */
+ integrator.type = &smc_double; /* Type of the Monte Carlo weight */
+ integrator.max_steps = max_steps; /* Realization count */
+ /* Solve */
time_current(&begin);
- /* solve */
- SMC(device_create(NULL, NULL, 8, NULL, &smc));
SMC(solve(smc, &integrator, &ctx, &estimator));
time_current(&end);
- time_sub(&elapsed, &end, &begin);
- time_dump(&elapsed, TIME_MIN|TIME_SEC|TIME_MSEC, NULL, buf, sizeof(buf));
- S3D(scene_end_session(scene));
-
- if(ctx.exit_failure) {
- printf("Too many failures. The scene might not match the prerequisites: it must be\n"
- "closed and its normals must point *into* the volume.\n");
+ if(ctx.error) {
+ fprintf(stderr,
+ "Too many failures. The scene might not match the prerequisites:\n"
+ "it must be closed and its normals must point *into* the volume.\n");
goto error;
}
- /* get simulation results */
+ time_sub(&elapsed, &end, &begin);
+ time_dump(&elapsed, TIME_MIN|TIME_SEC|TIME_MSEC, NULL, buf, sizeof(buf));
+
+ /* Print the simulation results */
SMC(estimator_get_status(estimator, &estimator_status));
length = SMC_DOUBLE(estimator_status.E);
sig_length = SMC_DOUBLE(estimator_status.SE);
-
- printf("\n4V/S = %g ~ %g +/- %g\n# failures: %d/%lu\nElapsed time: %s\n",
- expected, length, sig_length, ctx.nfailures, max_steps, buf);
+ printf("4V/S = %g ~ %g +/- %g\n# failures: %d/%lu\nElapsed time: %s\n",
+ reference, length, sig_length, ctx.nfailures, max_steps, buf);
exit:
- /* clean data */
- SMC(device_ref_put(smc));
- SMC(estimator_ref_put(estimator));
- if(scene) {
- int mask;
- S3D(scene_get_session_mask(scene, &mask));
- if(mask) S3D(scene_end_session(scene));
- S3D(scene_ref_put(scene));
- }
-
+ /* Clean-up data */
+ if(smc) SMC(device_ref_put(smc));
+ if(estimator) SMC(estimator_ref_put(estimator));
+ S3D(scene_get_session_mask(scene, &mask));
+ if(mask) S3D(scene_end_session(scene));
return res;
+
error:
- res = RES_UNKNOWN_ERR;
goto exit;
}
-/* import shapes in a scene */
+/* Create a S3D scene from an obj in a scene */
static res_T
-import_obj(const char *filename, struct s3d_scene **scene) {
+import_obj(const char* filename, struct s3d_scene** out_scene)
+{
struct s3d_device* s3d = NULL;
+ struct s3d_scene* scene = NULL;
struct s3daw* s3daw = NULL;
size_t i, count;
- res_T res;
- const int VERBOSE = 1;
-
- res = s3d_device_create(NULL, NULL, VERBOSE, &s3d);
- if (res != RES_OK) goto end;
- res = s3daw_create(NULL, NULL, NULL, NULL, s3d, VERBOSE, &s3daw);
- if (res != RES_OK) goto end;
- res = s3daw_load(s3daw, filename);
- if (res != RES_OK) goto end;
- res = s3daw_get_shapes_count(s3daw, &count);
- if (res != RES_OK) goto end;
- res = s3d_scene_create(s3d, scene);
- if (res != RES_OK) goto end;
- for (i = 0; i < count; i ++) {
+ const int VERBOSE = 0;
+ res_T res = RES_OK;
+ ASSERT(out_scene);
+
+ #define CALL(Func) if(RES_OK != (res = Func)) goto error
+ CALL(s3d_device_create(NULL, NULL, VERBOSE, &s3d));
+ CALL(s3daw_create(NULL, NULL, NULL, NULL, s3d, VERBOSE, &s3daw));
+ CALL(s3daw_load(s3daw, filename));
+ CALL(s3daw_get_shapes_count(s3daw, &count));
+ CALL(s3d_scene_create(s3d, &scene));
+
+ FOR_EACH(i, 0, count) {
struct s3d_shape* shape;
- if (res != RES_OK) goto end;
- res = s3daw_get_shape(s3daw, i, &shape);
- if (res != RES_OK) goto end;
- res = s3d_scene_attach_shape(*scene, shape);
- if (res != RES_OK) goto end;
+ CALL(s3daw_get_shape(s3daw, i, &shape));
+ CALL(s3d_scene_attach_shape(scene, shape));
}
+ #undef CALL
- /* release memory */
- end:
- if (s3daw) res = s3daw_ref_put(s3daw);
- if (s3d) res = s3d_device_ref_put(s3d);
-
+exit:
+ /* Release memory */
+ if(s3daw) S3DAW(ref_put(s3daw));
+ if(s3d) S3D(device_ref_put(s3d));
+ *out_scene = scene;
return res;
+error:
+ if(scene) {
+ S3D(scene_ref_put(scene));
+ scene = NULL;
+ }
+ goto exit;
}
-int main(int argc, char *argv[]) {
- struct s3d_scene *scene;
- size_t nsteps = 0;
- double ks = 0;
+int
+main(int argc, char* argv[])
+{
+ struct s3d_scene* scene = NULL;
+ unsigned long nsteps = 10000;
+ double ks = 1.0;
res_T res = RES_OK;
+ int err = 0;
- /* check command arguments */
+ /* Check command arguments */
if(argc < 2 || argc > 4) {
- if(argc < 2)
- printf("Argument required\n");
- else
- printf("Too many arguments\n");
- printf("usage: %s OBJ_FILE [SAMPLES_COUNT [K_SCATTERING]]\n", argv[0]);
- return RES_BAD_ARG;
+ printf("Usage: %s OBJ_FILE [SAMPLES_COUNT [K_SCATTERING]]\n", argv[0]);
+ goto error;
}
- /* import file's content in the scene */
+ /* Import file's content in the scene */
res = import_obj(argv[1], &scene);
if(res != RES_OK) {
- printf("Invalid file '%s'\n", argv[1]);
- return res;
+ fprintf(stderr, "Couldn't import `%s'\n", argv[1]);
+ goto error;
}
- /* set number of realizations */
- if(argc < 3) {
- nsteps = 10000; /* default value */
- }
- else {
- long ns = atol(argv[2]);
- if(ns <= 0) {
- printf("Invalid number of steps `%s'\n", argv[2]);
- return RES_BAD_ARG;
+ /* Set number of realizations */
+ if(argc >= 3) {
+ res = cstr_to_ulong(argv[2], &nsteps);
+ if(res != RES_OK) {
+ fprintf(stderr, "Invalid number of steps `%s'\n", argv[2]);
+ goto error;
}
- nsteps = (size_t)ns;
}
- /* set ks */
- if(argc < 4) {
- ks = 1; /* default value */
- }
- else {
- ks = atof(argv[3]);
- if(ks <= 0) {
- printf("Invalid k-scattering value `%s'\n", argv[3]);
- return RES_BAD_ARG;
+ /* Set ks */
+ if(argc >= 4) {
+ res = cstr_to_double(argv[3], &ks);
+ if(res != RES_OK || ks <= 0) {
+ fprintf(stderr, "Invalid k-scattering value `%s'\n", argv[3]);
+ goto error;
}
}
- /* compute */
- CHECK(compute_4v_s(scene, nsteps, ks), RES_OK);
-
- /* check memory leaks */
- CHECK(mem_allocated_size(), 0);
+ res = compute_4v_s(scene, nsteps, ks);
+ if(res != RES_OK) {
+ fprintf(stderr, "Error in 4V/S integration\n");
+ goto error;
+ }
- return 0;
+exit:
+ if(scene) S3D(scene_ref_put(scene));
+ if(mem_allocated_size() != 0) {
+ fprintf(stderr, "Memory leaks: %lu Bytes\n",
+ (unsigned long)mem_allocated_size());
+ err = 1;
+ }
+ return err;
+error:
+ err = 1;
+ goto exit;
}
+
