commit 8d01b745cb962d05578f3dd50438747fdef98e51
parent 2138be840bccc66722a6650914a1a959e9b7ae34
Author: Vincent Forest <vincent.forest@meso-star.com>
Date: Tue, 3 May 2016 10:01:59 +0200
Clean-up and refactoring
Diffstat:
3 files changed, 43 insertions(+), 55 deletions(-)
diff --git a/src/s4vs.c b/src/s4vs.c
@@ -26,7 +26,6 @@
* The fact that you are presently reading this means that you have had
* 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>
@@ -40,16 +39,12 @@
static res_T
compute_4v_s(struct s3d_scene* scene, const size_t max_steps, const double ks)
{
- char buf[512];
- struct time begin, end, elapsed;
struct s4vs_context ctx;
- 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;
+ float S, V, reference;
int mask;
res_T res = RES_OK;
ASSERT(scene && max_steps > 0 && ks > 0);
@@ -76,6 +71,7 @@ compute_4v_s(struct s3d_scene* scene, const size_t max_steps, const double ks)
/* Initialize context for MC computation */
ctx.scene = scene;
ctx.ks = ks;
+ ctx.g = PI/4.0;
/* Setup Star-MC */
SMC(device_create(NULL, NULL, SMC_NTHREADS_DEFAULT, NULL, &smc));
@@ -85,11 +81,7 @@ compute_4v_s(struct s3d_scene* scene, const size_t max_steps, const double ks)
integrator.max_failures = max_steps / 1000; /* cancel if 0.1% of the realization fail */
/* Solve */
- time_current(&begin);
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));
/* Print the simulation results */
SMC(estimator_get_status(estimator, &estimator_status));
@@ -101,11 +93,12 @@ compute_4v_s(struct s3d_scene* scene, const size_t max_steps, const double ks)
(unsigned long)estimator_status.NF);
goto error;
}
- length = SMC_DOUBLE(estimator_status.E);
- sig_length = SMC_DOUBLE(estimator_status.SE);
- printf("4V/S = %g ~ %g +/- %g\n# failures: %lu/%lu\nElapsed time: %s\n",
- reference, length, sig_length,
- (unsigned long)estimator_status.NF, (unsigned long)max_steps, buf);
+ printf("4V/S = %g ~ %g +/- %g\n# failures: %lu/%lu\n",
+ reference,
+ SMC_DOUBLE(estimator_status.E),
+ SMC_DOUBLE(estimator_status.SE),
+ (unsigned long)estimator_status.NF,
+ (unsigned long)max_steps);
exit:
/* Clean-up data */
diff --git a/src/s4vs_realization.c b/src/s4vs_realization.c
@@ -40,70 +40,64 @@ s4vs_realization(void* out_length, struct ssp_rng* rng, void* context)
struct s4vs_context* ctx = (struct s4vs_context*)context;
struct s3d_attrib attrib;
struct s3d_primitive prim;
- float normal[3], dir[4], pos[3], range[2], st[2];
+ float normal[3], u[4], x[3], range[2], st[2];
struct s3d_hit hit;
- double length = 0;
- double lambda = 0;
- float u, v, w;
-
- /* Sample the scene surface to define the random walk starting point */
+ double w = 0;
+ double sigma = 0;
+ int keep_running = 0;
+ float r0, r1, r2;
/* 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));
+ r0 = ssp_rng_canonical_float(rng);
+ r1 = ssp_rng_canonical_float(rng);
+ r2 = ssp_rng_canonical_float(rng);
+ S3D(scene_sample(ctx->scene, r0, r1, r2, &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(pos, attrib.value);
+ f3_set(x, attrib.value);
/* Cosine weighted sampling of the hemisphere around the sampled normal */
- ssp_ran_hemisphere_cos(rng, normal, dir);
+ ssp_ran_hemisphere_cos(rng, normal, u);
/* Find the 1st hit from the sampled location along the sampled direction */
range[0] = 1.e-6f;
- range[1] = FLT_MAX;
+ range[1] = (float)INF;
do {
- S3D(scene_trace_ray(ctx->scene, pos, dir, range, &hit));
+ S3D(scene_trace_ray(ctx->scene, x, u, range, &hit));
range[0] += 1.e-6f;
} while(S3D_PRIMITIVE_EQ(&hit.prim, &prim)); /* Handle self-hit */
/* No intersection <=> numerical imprecision or geometry leakage */
- if(S3D_HIT_NONE(&hit))
- return RES_UNKNOWN_ERR; /* failed realization */
+ if(S3D_HIT_NONE(&hit)) return RES_UNKNOWN_ERR;
- /* Here the starting point and a propagation direction have been sampled */
- /* and we have determined the distance of the geometry in this direction */
+ keep_running = 1;
+ while(keep_running) { /* Here we go for the diffuse random walk */
+ sigma = ssp_ran_exp(rng, ctx->ks); /* Sample a length according to 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(sigma < hit.distance) {
+ int i;
+ FOR_EACH(i, 0, 3) x[i] = x[i] + (float)sigma*u[i];
+ ssp_ran_sphere_hg(rng, u, ctx->g, u);
- if(lambda >= hit.distance) {
- /* Lambda exceeds the geometry distance: stop the random walk */
- length += hit.distance;
- break;
- }
+ /* sample a new direction */
+ range[0] = 0;
+ range[1] = (float)INF;
+ S3D(scene_trace_ray(ctx->scene, x, u, range, &hit));
- /* 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));
+ w = w + sigma;
- /* 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 */
- if(S3D_HIT_NONE(&hit))
- return RES_UNKNOWN_ERR; /* failed realization */
+ /* No intersection <=> numerical imprecision or geometry leakage */
+ if(S3D_HIT_NONE(&hit)) return RES_UNKNOWN_ERR;
+
+ } else { /* Stop the random walk */
+ w = w + hit.distance;
+ keep_running = 0;
+ }
}
- SMC_DOUBLE(out_length) = length;
+ SMC_DOUBLE(out_length) = w;
return RES_OK;
}
diff --git a/src/s4vs_realization.h b/src/s4vs_realization.h
@@ -38,6 +38,7 @@ struct s4vs_context {
struct s3d_sampler* sampler;
struct s3d_scene* scene;
double ks;
+ double g;
};
/*******************************************************************************
