commit 5d3819e5511d23327c30ef11eda00736899e3443
parent 6559593a00cf159aa3551f3860cff332a9d1b981
Author: Christophe Coustet <christophe.coustet@meso-star.com>
Date: Fri, 24 Mar 2017 13:14:49 +0100
Merge remote-tracking branch 'origin/develop' into feature_outputs
Diffstat:
21 files changed, 986 insertions(+), 65 deletions(-)
diff --git a/cmake/CMakeLists.txt b/cmake/CMakeLists.txt
@@ -169,6 +169,8 @@ if(NOT NO_TEST)
new_test(test_ssol_solver5)
new_test(test_ssol_solver6)
new_test(test_ssol_solver7)
+ new_test(test_ssol_solver8)
+ new_test(test_ssol_solver9)
new_test(test_ssol_sun)
build_test(test_ssol_draw)
diff --git a/src/ssol.h b/src/ssol.h
@@ -66,6 +66,12 @@ enum ssol_side_flag {
SSOL_INVALID_SIDE = BIT(2)
};
+enum ssol_path_type {
+ SSOL_PATH_MISSING, /* The path misses the receivers */
+ SSOL_PATH_SHADOW, /* The path is occluded before the sampled geometry */
+ SSOL_PATH_SUCCESS /* The path contributes to at least one receiver */
+};
+
enum ssol_material_type {
SSOL_MATERIAL_MATTE,
SSOL_MATERIAL_MIRROR,
@@ -298,6 +304,27 @@ struct ssol_instantiated_shaded_shape {
static const struct ssol_instantiated_shaded_shape
SSOL_INSTANTIATED_SHADED_SHAPE_NULL = SSOL_INSTANTIATED_SHADED_SHAPE_NULL__;
+struct ssol_path_tracker {
+ /* Control the length of the path segment starting/ending from/to the
+ * infinite. A value less than zero means for default value */
+ double sun_ray_length;
+ double infinite_ray_length;
+};
+
+#define SSOL_PATH_TRACKER_DEFAULT__ {-1, -1}
+static const struct ssol_path_tracker SSOL_PATH_TRACKER_DEFAULT =
+ SSOL_PATH_TRACKER_DEFAULT__;
+
+struct ssol_path {
+ /* Internal data */
+ const void* path__;
+};
+
+struct ssol_path_vertex {
+ double pos[3]; /* Position */
+ double weight; /* Monte-Carlo weight */
+};
+
struct ssol_mc_result {
double E; /* Expectation */
double V; /* Variance */
@@ -998,6 +1025,36 @@ ssol_estimator_get_mc_sampled
struct ssol_mc_sampled* sampled);
/*******************************************************************************
+ * Tracked paths
+ ******************************************************************************/
+SSOL_API res_T
+ssol_estimator_get_tracked_paths_count
+ (const struct ssol_estimator* estimator,
+ size_t* npaths);
+
+SSOL_API res_T
+ssol_estimator_get_tracked_path
+ (const struct ssol_estimator* estimator,
+ const size_t ipath,
+ struct ssol_path* path);
+
+SSOL_API res_T
+ssol_path_get_vertices_count
+ (const struct ssol_path* path,
+ size_t* nvertices);
+
+SSOL_API res_T
+ssol_path_get_vertex
+ (const struct ssol_path* path,
+ const size_t ivertex,
+ struct ssol_path_vertex* vertex);
+
+SSOL_API res_T
+ssol_path_get_type
+ (const struct ssol_path* path,
+ enum ssol_path_type* type);
+
+/*******************************************************************************
* Per receiver MC estimations
******************************************************************************/
SSOL_API res_T
@@ -1027,6 +1084,7 @@ ssol_solve
(struct ssol_scene* scn,
struct ssp_rng* rng,
const size_t realisations_count,
+ const struct ssol_path_tracker* tracker, /* NULL<=>Do not record the paths */
FILE* output, /* May be NULL <=> does not ouput ssol_receiver_data */
struct ssol_estimator** estimator);
diff --git a/src/ssol_draw_pt.c b/src/ssol_draw_pt.c
@@ -104,7 +104,12 @@ sun_lighting
ASSERT(sun && view && ray_data && bsdf && wo && N && ray_org);
ASSERT(d3_dot(wo, N) >= 0); /* Assume that wo point outward the surface */
+ /* Ensure that the incoming direction point outward the surface */
d3_minus(wi, sun->direction);
+
+ /* The point look backward the sun */
+ if(d3_dot(wi, N) < 0) return 0.0;
+
R = ssf_bsdf_eval(bsdf, wo, N, wi);
if(R <= 0) return 0.0;
diff --git a/src/ssol_estimator.c b/src/ssol_estimator.c
@@ -20,9 +20,10 @@
#include "ssol_device_c.h"
#include "ssol_instance_c.h"
-#include <rsys/rsys.h>
+#include <rsys/double3.h>
#include <rsys/mem_allocator.h>
#include <rsys/ref_count.h>
+#include <rsys/rsys.h>
#include <math.h>
@@ -79,6 +80,7 @@ estimator_release(ref_T* ref)
dev = estimator->dev;
htable_receiver_release(&estimator->mc_receivers);
htable_sampled_release(&estimator->mc_sampled);
+ darray_path_release(&estimator->paths);
ASSERT(dev && dev->allocator);
MEM_RM(dev->allocator, estimator);
SSOL(device_ref_put(dev));
@@ -233,6 +235,60 @@ ssol_estimator_get_mc_sampled
return RES_OK;
}
+res_T
+ssol_estimator_get_tracked_paths_count
+ (const struct ssol_estimator* estimator, size_t* npaths)
+{
+ if(!estimator || !npaths) return RES_BAD_ARG;
+ *npaths = darray_path_size_get(&estimator->paths);
+ return RES_OK;
+}
+
+res_T
+ssol_estimator_get_tracked_path
+ (const struct ssol_estimator* estimator,
+ const size_t ipath,
+ struct ssol_path* path)
+{
+ if(!estimator || ipath >= darray_path_size_get(&estimator->paths) || !path)
+ return RES_BAD_ARG;
+ path->path__ = darray_path_cdata_get(&estimator->paths) + ipath;
+ return RES_OK;
+}
+
+res_T
+ssol_path_get_vertices_count(const struct ssol_path* path, size_t* nvertices)
+{
+ const struct path* p;
+ if(!path || !nvertices) return RES_BAD_ARG;
+ p = path->path__;
+ *nvertices = darray_path_vertex_size_get(&p->vertices);
+ return RES_OK;
+}
+
+res_T
+ssol_path_get_vertex
+ (const struct ssol_path* path,
+ const size_t ivertex,
+ struct ssol_path_vertex* vertex)
+{
+ const struct path* p;
+ if(!path || !vertex) return RES_BAD_ARG;
+ p = path->path__;
+ if(ivertex >= darray_path_vertex_size_get(&p->vertices)) return RES_BAD_ARG;
+ *vertex = darray_path_vertex_cdata_get(&p->vertices)[ivertex];
+ return RES_OK;
+}
+
+res_T
+ssol_path_get_type(const struct ssol_path* path, enum ssol_path_type* type)
+{
+ ASSERT(path && type);
+ if(!path || !type) return RES_BAD_ARG;
+ *type = ((struct path*)path->path__)->type;
+ return RES_OK;
+}
+
/*******************************************************************************
* Local function
******************************************************************************/
@@ -258,6 +314,7 @@ estimator_create
htable_receiver_init(dev->allocator, &estimator->mc_receivers);
htable_sampled_init(dev->allocator, &estimator->mc_sampled);
+ darray_path_init(dev->allocator, &estimator->paths);
SSOL(device_ref_get(dev));
estimator->dev = dev;
ref_init(&estimator->ref);
diff --git a/src/ssol_estimator_c.h b/src/ssol_estimator_c.h
@@ -365,6 +365,77 @@ error:
#include <rsys/hash_table.h>
/*******************************************************************************
+ * Radiative path
+ ******************************************************************************/
+#define DARRAY_NAME path_vertex
+#define DARRAY_DATA struct ssol_path_vertex
+#include <rsys/dynamic_array.h>
+
+struct path {
+ enum ssol_path_type type;
+ struct darray_path_vertex vertices;
+};
+
+static INLINE void
+path_init(struct mem_allocator* allocator, struct path* path)
+{
+ ASSERT(path);
+ path->type = SSOL_PATH_MISSING;
+ darray_path_vertex_init(allocator, &path->vertices);
+}
+
+static INLINE void
+path_release(struct path* path)
+{
+ ASSERT(path);
+ darray_path_vertex_release(&path->vertices);
+}
+
+static INLINE res_T
+path_copy(struct path* dst, const struct path* src)
+{
+ ASSERT(dst && src);
+ dst->type = src->type;
+ return darray_path_vertex_copy(&dst->vertices, &src->vertices);
+}
+
+static INLINE res_T
+path_copy_and_release(struct path* dst, struct path* src)
+{
+ ASSERT(dst && src);
+ dst->type = src->type;
+ return darray_path_vertex_copy_and_release(&dst->vertices, &src->vertices);
+}
+
+static INLINE res_T
+path_copy_and_clear(struct path* dst, struct path* src)
+{
+ ASSERT(dst && src);
+ dst->type = src->type;
+ return darray_path_vertex_copy_and_clear(&dst->vertices, &src->vertices);
+}
+
+static INLINE res_T
+path_add_vertex(struct path* path, const double pos[3], const double weight)
+{
+ struct ssol_path_vertex vertex;
+ ASSERT(path && pos && weight >= 0);
+ vertex.pos[0] = pos[0];
+ vertex.pos[1] = pos[1];
+ vertex.pos[2] = pos[2];
+ vertex.weight = weight;
+ return darray_path_vertex_push_back(&path->vertices, &vertex);
+}
+
+#define DARRAY_NAME path
+#define DARRAY_DATA struct path
+#define DARRAY_FUNCTOR_INIT path_init
+#define DARRAY_FUNCTOR_RELEASE path_release
+#define DARRAY_FUNCTOR_COPY path_copy
+#define DARRAY_FUNCTOR_COPY_AND_RELEASE path_copy_and_release
+#include <rsys/dynamic_array.h>
+
+/*******************************************************************************
* Estimator data structure
******************************************************************************/
struct ssol_estimator {
@@ -379,6 +450,8 @@ struct ssol_estimator {
struct htable_receiver mc_receivers; /* Per receiver MC */
struct htable_sampled mc_sampled; /* Per sampled instance MC */
+ struct darray_path paths; /* Tracked paths */
+
/* Overall area of the sampled instances. Actually this is not the area that
* is effectively sampled since an instance may be sampled through a proxy
* geometry */
diff --git a/src/ssol_solver.c b/src/ssol_solver.c
@@ -55,6 +55,7 @@ struct thread_context {
struct mc_data cos_factor;
struct htable_receiver mc_rcvs;
struct htable_sampled mc_samps;
+ struct darray_path paths; /* paths */
};
static void
@@ -65,6 +66,7 @@ thread_context_release(struct thread_context* ctx)
if(ctx->bsdf) SSF(bsdf_ref_put(ctx->bsdf));
htable_receiver_release(&ctx->mc_rcvs);
htable_sampled_release(&ctx->mc_samps);
+ darray_path_release(&ctx->paths);
}
static res_T
@@ -76,6 +78,7 @@ thread_context_init(struct mem_allocator* allocator, struct thread_context* ctx)
memset(ctx, 0, sizeof(ctx[0]));
htable_receiver_init(allocator, &ctx->mc_rcvs);
htable_sampled_init(allocator, &ctx->mc_samps);
+ darray_path_init(allocator, &ctx->paths);
res = ssf_bsdf_create(allocator, &ctx->bsdf);
if(res != RES_OK) goto error;
@@ -104,6 +107,8 @@ thread_context_copy
if(res != RES_OK) return res;
res = htable_sampled_copy(&dst->mc_samps, &src->mc_samps);
if(res != RES_OK) return res;
+ res = darray_path_copy(&dst->paths, &src->paths);
+ if(res != RES_OK) return res;
return RES_OK;
}
@@ -114,6 +119,7 @@ thread_context_clear(struct thread_context* ctx)
if(ctx->rng) SSP(rng_ref_put(ctx->rng));
htable_receiver_clear(&ctx->mc_rcvs);
htable_sampled_clear(&ctx->mc_samps);
+ darray_path_clear(&ctx->paths);
}
static res_T
@@ -196,7 +202,6 @@ point_init
struct s3d_hit hit;
struct ray_data ray_data = RAY_DATA_NULL;
float dir[3], pos[3], range[2] = { 0, FLT_MAX };
- double cos_sun;
size_t id;
res_T res = RES_OK;
ASSERT(pt && scn && sampled && view_samp && view_rt);
@@ -219,15 +224,7 @@ point_init
f3_normalize(attr.value, attr.value);
d3_set_f3(pt->N, attr.value);
- /* Sample a sun direction */
- ranst_sun_dir_get(ran_sun_dir, rng, pt->dir);
-
- /* Initialise the Monte Carlo weight */
- cos_sun = fabs(d3_dot(pt->N, pt->dir));
- pt->weight = scn->sun->dni * sampled_area_proxy * cos_sun;
- pt->absorptivity_loss = pt->reflectivity_loss = 0;
pt->cos_factor = cos_sun;
-
/* Retrieve the sampled instance and shaded shape */
pt->inst = *htable_instance_find(&scn->instances_samp, &pt->prim.inst_id);
id = *htable_shaded_shape_find
@@ -235,17 +232,34 @@ point_init
pt->sshape = darray_shaded_shape_cdata_get
(&pt->inst->object->shaded_shapes) + id;
- /* Store sampled entity related weights */
- res = get_mc_sampled(sampled, pt->inst, &pt->mc_samp);
- if(res != RES_OK) goto error;
- pt->mc_samp->nb_samples++;
+ /* Sample a sun direction */
+ ranst_sun_dir_get(ran_sun_dir, rng, pt->dir);
- /* For punched surface, retrieve the sampled position and normal onto the
- * quadric surface */
+ /* Initialise the Monte Carlo weight */
if(pt->sshape->shape->type == SHAPE_PUNCHED) {
+ double proxy_sun_cos = fabs(d3_dot(pt->N, pt->dir));
+ double cos_ratio, surface_proxy_cos, surface_sun_cos, tmp_n[3];
+ /* For punched surface, retrieve the sampled position and normal onto the
+ * quadric surface */
punched_shape_project_point
- (pt->sshape->shape, pt->inst->transform, pt->pos, pt->pos, pt->N);
+ (pt->sshape->shape, pt->inst->transform, pt->pos, pt->pos, tmp_n);
+ surface_proxy_cos = d3_dot(pt->N, tmp_n);
+ surface_sun_cos = d3_dot(tmp_n, pt->dir);
+ cos_ratio = fabs(surface_sun_cos / surface_proxy_cos);
+ d3_set(pt->N, tmp_n);
+ pt->weight = scn->sun->dni * sampled_area_proxy * cos_ratio;
+ pt->cos_loss = scn->sun->dni * sampled_area_proxy * (1 - proxy_sun_cos);
+ } else {
+ double surface_sun_cos = fabs(d3_dot(pt->N, pt->dir));
+ pt->weight = scn->sun->dni * sampled_area_proxy * surface_sun_cos;
+ pt->cos_loss = scn->sun->dni * sampled_area_proxy * (1 - surface_sun_cos);
}
+ pt->absorptivity_loss = pt->reflectivity_loss = 0;
+
+ /* Store sampled entity related weights */
+ res = get_mc_sampled(sampled, pt->inst, &pt->mc_samp);
+ if(res != RES_OK) goto error;
+ pt->mc_samp->nb_samples++;
/* Define the primitive side on which the point lies */
if(d3_dot(pt->N, pt->dir) < 0) {
@@ -460,6 +474,36 @@ check_scene(const struct ssol_scene* scene, const char* caller)
return RES_OK;
}
+/* Compute an empirical length of the path segment coming from/going to the
+ * infinite, wrt the scene bounding box */
+static INLINE double
+compute_infinite_path_segment_extend(struct s3d_scene_view* view)
+{
+ float lower[3], upper[3], size[3];
+ ASSERT(view);
+ S3D(scene_view_get_aabb(view, lower, upper));
+ f3_sub(size, upper, lower);
+ return MMAX(size[0], MMAX(size[1], size[2])) * 0.75;
+}
+
+static INLINE res_T
+path_register_and_clear
+ (struct darray_path* paths,
+ struct path* path)
+{
+ struct path* dst_path;
+ size_t ipath;
+ res_T res = RES_OK;
+ ASSERT(paths && path);
+
+ ipath = darray_path_size_get(paths);
+ res = darray_path_resize(paths, ipath + 1);
+ if(res != RES_OK) return res;
+
+ dst_path = darray_path_data_get(paths) + ipath;
+ return path_copy_and_clear(dst_path, path);
+}
+
static res_T
accum_mc_receivers_1side
(struct mc_receiver_1side* dst,
@@ -632,7 +676,7 @@ update_mc
res = mc_shape_1side_get_mc_primitive(mc_shape1, pt->prim.prim_id, &mc_prim1);
if(res != RES_OK) goto error;
- #define ACCUM_WEIGHT(Name, W) { \
+ #define ACCUM_WEIGHT(Name, W) { \
mc_prim1->Name.weight += (W); \
mc_prim1->Name.sqr_weight += (W)*(W); \
} (void)0
@@ -659,8 +703,10 @@ trace_radiative_path
struct s3d_scene_view* view_rt,
struct ranst_sun_dir* ran_sun_dir,
struct ranst_sun_wl* ran_sun_wl,
+ const struct ssol_path_tracker* tracker, /* May be NULL */
FILE* output) /* May be NULL */
{
+ struct path path;
struct s3d_hit hit = S3D_HIT_NULL;
struct point pt;
float org[3], dir[3], range[2] = { 0, FLT_MAX };
@@ -669,20 +715,40 @@ trace_radiative_path
res_T res = RES_OK;
ASSERT(thread_ctx && scn && view_samp && view_rt && ran_sun_dir && ran_sun_wl);
+ if(tracker) path_init(scn->dev->allocator, &path);
+
/* Find a new starting point of the radiative random walk */
res = point_init(&pt, sampled_area_proxy, scn, &thread_ctx->mc_samps,
view_samp, view_rt, ran_sun_dir, ran_sun_wl, thread_ctx->rng, &is_lit);
if(res != RES_OK) goto error;
+ if(tracker) {
+ /* Add the first point of the starting segment */
+ if(tracker->sun_ray_length > 0) {
+ double pos[3], wi[3];
+ d3_minus(wi, pt.dir);
+ d3_muld(wi, wi, tracker->sun_ray_length);
+ d3_add(pos, pt.pos, wi);
+ res = path_add_vertex(&path, pos, scn->sun->dni);
+ if(res != RES_OK) goto error;
+ }
+
+ /* Register the init position onto the sampled geometry */
+ res = path_add_vertex(&path, pt.pos, pt.weight);
+ if(res != RES_OK) goto error;
+ }
+
#define ACCUM_WEIGHT(Res, W) { \
Res.weight += (W); \
Res.sqr_weight += (W)*(W); \
} (void)0
ACCUM_WEIGHT(thread_ctx->cos_factor, pt.cos_factor);
+
if(!is_lit) { /* The starting point is not lit */
ACCUM_WEIGHT(pt.mc_samp->shadowed, pt.weight);
ACCUM_WEIGHT(thread_ctx->shadowed, pt.weight);
#undef ACCUM_WEIGHT
+ if(tracker) path.type = SSOL_PATH_SHADOW;
} else {
int hit_a_receiver = 0;
@@ -740,7 +806,17 @@ trace_radiative_path
ray_data.range_min = range[0];
ray_data.dst = FLT_MAX;
S3D(scene_view_trace_ray(view_rt, org, dir, range, &ray_data, &hit));
- if(S3D_HIT_NONE(&hit)) break;
+ if(S3D_HIT_NONE(&hit)) {
+ /* Add the point of the last path segment going to the infinite */
+ if(tracker && tracker->infinite_ray_length > 0) {
+ double pos[3], wi[3];
+ d3_set_f3(wi, dir);
+ d3_add(pos, pt.pos, d3_muld(wi, wi, tracker->infinite_ray_length));
+ res = path_add_vertex(&path, pos, pt.weight);
+ if(res != RES_OK) goto error;
+ }
+ break;
+ }
depth += mtl->type != SSOL_MATERIAL_VIRTUAL;
@@ -755,14 +831,28 @@ trace_radiative_path
/* Update the point */
point_update_from_hit(&pt, scn, org, dir, &hit, &ray_data);
+
+ if(tracker) {
+ res = path_add_vertex(&path, pt.pos, pt.weight);
+ if(res != RES_OK) goto error;
+ }
}
if(!hit_a_receiver) {
thread_ctx->missing.weight += pt.weight;
thread_ctx->missing.sqr_weight += pt.weight*pt.weight;
}
+ if(tracker) {
+ path.type = hit_a_receiver ? SSOL_PATH_SUCCESS : SSOL_PATH_MISSING;
+ }
+ }
+
+ if(tracker) {
+ res = path_register_and_clear(&thread_ctx->paths, &path);
+ if(res != RES_OK) goto error;
}
exit:
+ if(tracker) path_release(&path);
return res;
error:
goto exit;
@@ -776,6 +866,7 @@ ssol_solve
(struct ssol_scene* scn,
struct ssp_rng* rng_state,
const size_t realisations_count,
+ const struct ssol_path_tracker* path_tracker,
FILE* output,
struct ssol_estimator** out_estimator)
{
@@ -787,6 +878,7 @@ ssol_solve
struct ranst_sun_wl* ran_sun_wl = NULL;
struct darray_thread_ctx thread_ctxs;
struct ssol_estimator* estimator = NULL;
+ struct ssol_path_tracker tracker;
struct ssp_rng_proxy* rng_proxy = NULL;
double sampled_area;
double sampled_area_proxy;
@@ -839,6 +931,17 @@ ssol_solve
if(res != RES_OK) goto error;
}
+ /* Setup the path tracker */
+ if(path_tracker) {
+ tracker = *path_tracker;
+ if(tracker.sun_ray_length < 0 || tracker.infinite_ray_length < 0) {
+ const double extend = compute_infinite_path_segment_extend(view_rt);
+ if(tracker.sun_ray_length < 0) tracker.sun_ray_length = extend;
+ if(tracker.infinite_ray_length < 0) tracker.infinite_ray_length = extend;
+ }
+ path_tracker = &tracker;
+ }
+
/* Launch the parallel MC estimation */
#pragma omp parallel for schedule(static)
for(i = 0; i < nrealisations; ++i) {
@@ -853,7 +956,7 @@ ssol_solve
/* Execute a MC experiment */
res_local = trace_radiative_path((size_t)i, sampled_area_proxy, thread_ctx,
- scn, view_samp, view_rt, ran_sun_dir, ran_sun_wl, output);
+ scn, view_samp, view_rt, ran_sun_dir, ran_sun_wl, path_tracker, output);
if(res_local != RES_OK) {
ATOMIC_SET(&res, res_local);
continue;
@@ -922,6 +1025,23 @@ ssol_solve
}
}
+ /* Merge per thread tracked paths */
+ if(path_tracker) {
+ FOR_EACH(i, 0, nthreads) {
+ struct thread_context* thread_ctx;
+ size_t ipath, npaths;
+
+ thread_ctx = darray_thread_ctx_data_get(&thread_ctxs) + i;
+ npaths = darray_path_size_get(&thread_ctx->paths);
+ FOR_EACH(ipath, 0, npaths) {
+ struct path* path;
+ path = darray_path_data_get(&thread_ctx->paths) + ipath;
+ res = path_register_and_clear(&estimator->paths, path);
+ if(res != RES_OK) goto error;
+ }
+ }
+ }
+
estimator->realisation_count = realisations_count;
estimator->sampled_area = sampled_area;
diff --git a/src/test_ssol_by_receiver_integration.c b/src/test_ssol_by_receiver_integration.c
@@ -134,19 +134,19 @@ main(int argc, char** argv)
#define S_DNI_cos (4 * 1000 * cos(PI / 4))
#define GET_MC_RCV ssol_estimator_get_mc_receiver
#define GET_MC_SAMP_X_RCV ssol_estimator_get_mc_sampled_x_receiver
- CHECK(ssol_solve(scene, rng, N__, NULL, &estimator1), RES_OK);
+ CHECK(ssol_solve(scene, rng, N__, 0, NULL, &estimator1), RES_OK);
CHECK(GET_MC_RCV(estimator1, target, SSOL_FRONT, &mc_rcv), RES_OK);
printf("Ir(target) = %g +/- %g\n",
mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE);
CHECK(ssol_instance_set_receiver(heliostat, SSOL_FRONT, 0), RES_OK);
CHECK(eq_eps(mc_rcv.integrated_irradiance.E, S_DNI_cos, S_DNI_cos * 2e-1), 1);
- CHECK(ssol_solve(scene, rng, 8 * N__, NULL, &estimator2), RES_OK);
+ CHECK(ssol_solve(scene, rng, 8 * N__, 0, NULL, &estimator2), RES_OK);
CHECK(GET_MC_RCV(estimator2, target, SSOL_FRONT, &mc_rcv), RES_OK);
printf("Ir(target) = %g +/- %g\n",
mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE);
CHECK(eq_eps(mc_rcv.integrated_irradiance.E, S_DNI_cos, S_DNI_cos * 5e-2), 1);
CHECK(ssol_estimator_ref_put(estimator1), RES_OK);
- CHECK(ssol_solve(scene, rng, 3 * N__, NULL, &estimator1), RES_OK);
+ CHECK(ssol_solve(scene, rng, 3 * N__, 0, NULL, &estimator1), RES_OK);
CHECK(GET_MC_RCV(estimator1, target, SSOL_FRONT, &mc_rcv), RES_OK);
printf("Ir(target) = %g +/- %g\n",
mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE);
diff --git a/src/test_ssol_circ2D_geometry.h b/src/test_ssol_circ2D_geometry.h
@@ -0,0 +1,54 @@
+/* Copyright (C) CNRS 2016-2017
+ *
+ * This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. */
+
+#include "test_ssol_geometries.h"
+#include <rsys/math.h>
+
+
+/*******************************************************************************
+* Circle polygon
+******************************************************************************/
+#if !defined(RADIUS)
+#error "Missing the RADIUS macro defining the circle radius"
+#endif
+#if !defined(NVERTS)
+#define NVERTS 36
+#endif
+#if !defined(CIRCLE_NAME)
+#error "Missing the CIRCLE_NAME macro defining the circle name"
+#endif
+
+#define EDGES__ CONCAT(CIRCLE_NAME, _EDGES__)
+#define INIT_FUNC__ CONCAT(init_, CIRCLE_NAME)
+
+/* should be const but scpr expects non-const data */
+static double EDGES__ [2*NVERTS];
+
+static void INIT_FUNC__() {
+ int n;
+ /* radius that give the same area than a circle */
+ double r = sqrt(2 * RADIUS * RADIUS * PI / (sin(2 * PI / NVERTS) * NVERTS));
+ for (n = 0; n < NVERTS; n++) {
+ EDGES__[2 * n] = r * cos((double)-n * 2 * PI / (double)NVERTS);
+ EDGES__[2 * n + 1] = r * sin((double)-n * 2 * PI / (double) NVERTS);
+ }
+}
+
+#undef EDGES__
+#undef INIT_FUNC__
+
+#undef RADIUS
+#undef NVERTS
+#undef CIRCLE_NAME
diff --git a/src/test_ssol_cube_geometry.h b/src/test_ssol_cube_geometry.h
@@ -0,0 +1,107 @@
+/* Copyright (C) CNRS 2016-2017
+ *
+ * This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. */
+
+#include "test_ssol_geometries.h"
+
+/*******************************************************************************
+* Rectangle polygon
+******************************************************************************/
+#if !defined(HALF_X) && !(defined(X_MIN) && defined(X_MAX))
+#error "Missing the HALF_X or X_MIN and X_MAX macros defining the cube size"
+#endif
+#if !defined(HALF_Y) && !(defined(Y_MIN) && defined(Y_MAX))
+#error "Missing the HALF_Y or Y_MIN and Y_MAX macros defining the cube size"
+#endif
+#if !defined(HALF_Z) && !(defined(Z_MIN) && defined(Z_MAX))
+#error "Missing the HALF_Z or Z_MIN and Z_MAX macros defining the cube size"
+#endif
+#if !defined(CUBE_NAME)
+#error "Missing the CUBE_NAME macro defining the rectangle name"
+#endif
+
+#define EDGES__ CONCAT(CUBE_NAME, _EDGES__)
+#define CUBE_NVERTS__ CONCAT(CUBE_NAME, _NVERTS__)
+
+#if !defined(X_MIN)
+#define X_MIN (float)(-(HALF_X))
+#endif
+
+#if !defined(X_MAX)
+#define X_MAX (float)(HALF_X)
+#endif
+
+#if !defined(Y_MIN)
+#define Y_MIN (float)(-(HALF_Y))
+#endif
+
+#if !defined(Y_MAX)
+#define Y_MAX (float)(HALF_Y)
+#endif
+
+#if !defined(Z_MIN)
+#define Z_MIN (float)(-(HALF_Z))
+#endif
+
+#if !defined(Z_MAX)
+#define Z_MAX (float)(HALF_Z)
+#endif
+
+static const float EDGES__ [] = {
+ X_MIN, Y_MIN, Z_MIN,
+ X_MIN, Y_MIN, Z_MAX,
+ X_MIN, Y_MAX, Z_MIN,
+ X_MIN, Y_MAX, Z_MAX,
+ X_MAX, Y_MIN, Z_MIN,
+ X_MAX, Y_MIN, Z_MAX,
+ X_MAX, Y_MAX, Z_MIN,
+ X_MAX, Y_MAX, Z_MAX
+};
+
+const unsigned CUBE_NVERTS__ = sizeof(EDGES__) / sizeof(float[3]);
+
+const unsigned TRG_IDS__ [] = {
+ 0, 6, 4,
+ 0, 2, 6,
+ 0, 3, 2,
+ 0, 1, 3,
+ 2, 7, 6,
+ 2, 3, 7,
+ 4, 6, 7,
+ 4, 7, 5,
+ 0, 4, 5,
+ 0, 5, 1,
+ 1, 5, 7,
+ 1, 7, 3
+};
+const unsigned CUBE_NTRIS__ = sizeof(TRG_IDS__) / sizeof(unsigned[3]);
+
+static const struct desc CUBE_DESC__ = { EDGES__, TRG_IDS__ };
+
+#undef EDGES__
+#undef TRG_IDS__
+#undef CUBE_DESC__
+#undef CUBE_NVERTS__
+#undef CUBE_NTRIS__
+
+#undef HALF_X
+#undef HALF_Y
+#undef HALF_Z
+#undef X_MIN
+#undef X_MAX
+#undef Y_MIN
+#undef Y_MAX
+#undef Z_MIN
+#undef Z_MAX
+#undef CUBE_NAME
diff --git a/src/test_ssol_rect2D_geometry.h b/src/test_ssol_rect2D_geometry.h
@@ -18,25 +18,41 @@
/*******************************************************************************
* Rectangle polygon
******************************************************************************/
-#if !defined(HALF_X)
-#error "Missing the HALF_X macro defining the rectangle size"
+#if !defined(HALF_X) && !(defined(X_MIN) && defined(X_MAX))
+#error "Missing the HALF_X or X_MIN and X_MAX macros defining the rectangle size"
#endif
-#if !defined(HALF_Y)
-#error "Missing the HALF_Y macro defining the rectangle size"
+#if !defined(HALF_Y) && !(defined(Y_MIN) && defined(Y_MAX))
+#error "Missing the HALF_Y or Y_MIN and Y_MAX macros defining the rectangle size"
#endif
#if !defined(POLYGON_NAME)
-#error "Missing the POLYGON_NAME macro defining the polygon name"
+#error "Missing the POLYGON_NAME macro defining the rectangle name"
#endif
#define EDGES__ CONCAT(POLYGON_NAME, _EDGES__)
#define POLY_NVERTS__ CONCAT(POLYGON_NAME, _NVERTS__)
+#if !defined(X_MIN)
+#define X_MIN (float)(-(HALF_X))
+#endif
+
+#if !defined(X_MAX)
+#define X_MAX (float)(HALF_X)
+#endif
+
+#if !defined(Y_MIN)
+#define Y_MIN (float)(-(HALF_Y))
+#endif
+
+#if !defined(Y_MAX)
+#define Y_MAX (float)(HALF_Y)
+#endif
+
/* should be const but scpr expects non-const data */
static double EDGES__ [] = {
- -HALF_X, -HALF_Y,
- -HALF_X, HALF_Y,
- HALF_X, HALF_Y,
- HALF_X, -HALF_Y,
+ X_MIN, Y_MIN,
+ X_MIN, Y_MAX,
+ X_MAX, Y_MAX,
+ X_MAX, Y_MIN
};
const unsigned POLY_NVERTS__ = sizeof(EDGES__) / sizeof(double[2]);
@@ -46,4 +62,8 @@ const unsigned POLY_NVERTS__ = sizeof(EDGES__) / sizeof(double[2]);
#undef HALF_X
#undef HALF_Y
+#undef X_MIN
+#undef X_MAX
+#undef Y_MIN
+#undef Y_MAX
#undef POLYGON_NAME
diff --git a/src/test_ssol_rect_geometry.h b/src/test_ssol_rect_geometry.h
@@ -18,11 +18,11 @@
/*******************************************************************************
* Rectangle plane
******************************************************************************/
-#if !defined(HALF_X)
-#error "Missing the HALF_X macro defining the rectangle size"
+#if !defined(HALF_X) && !(defined(X_MIN) && defined(X_MAX))
+#error "Missing the HALF_X or X_MIN and X_MAX macros defining the rectangle size"
#endif
-#if !defined(HALF_Y)
-#error "Missing the HALF_Y macro defining the rectangle size"
+#if !defined(HALF_Y) && !(defined(Y_MIN) && defined(Y_MAX))
+#error "Missing the HALF_Y or Y_MIN and Y_MAX macros defining the rectangle size"
#endif
#if !defined(PLANE_NAME)
#error "Missing the DARRAY_NAME macro defining the rectangle name"
@@ -34,11 +34,27 @@
#define RECT_NVERTS__ CONCAT(PLANE_NAME, _NVERTS__)
#define RECT_NTRIS__ CONCAT(PLANE_NAME, _NTRIS__)
+#if !defined(X_MIN)
+#define X_MIN (float)(-(HALF_X))
+#endif
+
+#if !defined(X_MAX)
+#define X_MAX (float)(HALF_X)
+#endif
+
+#if !defined(Y_MIN)
+#define Y_MIN (float)(-(HALF_Y))
+#endif
+
+#if !defined(Y_MAX)
+#define Y_MAX (float)(HALF_Y)
+#endif
+
static const float EDGES__ [] = {
- (float) -HALF_X, (float) -HALF_Y, 0.f,
- (float) HALF_X, (float) -HALF_Y, 0.f,
- (float) HALF_X, (float) HALF_Y, 0.f,
- (float) -HALF_X, (float) HALF_Y, 0.f
+ X_MIN, Y_MIN, 0.f,
+ X_MAX, Y_MIN, 0.f,
+ X_MAX, Y_MAX, 0.f,
+ X_MIN, Y_MAX, 0.f
};
const unsigned RECT_NVERTS__ = sizeof(EDGES__) / sizeof(float[3]);
@@ -56,4 +72,8 @@ static const struct desc RECT_DESC__ = { EDGES__, TRG_IDS__ };
#undef HALF_X
#undef HALF_Y
+#undef X_MIN
+#undef X_MAX
+#undef Y_MIN
+#undef Y_MAX
#undef PLANE_NAME
diff --git a/src/test_ssol_solver1.c b/src/test_ssol_solver1.c
@@ -74,6 +74,8 @@ main(int argc, char** argv)
struct ssol_mc_receiver mc_rcv;
struct ssol_mc_shape mc_shape;
struct ssol_mc_primitive mc_prim;
+ struct ssol_path path;
+ struct ssol_path_vertex vertex;
double dir[3];
double wavelengths[3] = { 1, 2, 3 };
double intensities[3] = { 1, 0.8, 1 };
@@ -121,14 +123,14 @@ main(int argc, char** argv)
CHECK(ssol_scene_create(dev, &scene), RES_OK);
CHECK(ssol_scene_attach_sun(scene, sun), RES_OK);
- CHECK(ssol_solve(NULL, rng, 10, NULL, &estimator), RES_BAD_ARG);
- CHECK(ssol_solve(scene, NULL, 10, NULL, &estimator), RES_BAD_ARG);
- CHECK(ssol_solve(scene, rng, 0, NULL, &estimator), RES_BAD_ARG);
- CHECK(ssol_solve(scene, rng, 10, NULL, &estimator), RES_BAD_ARG);
- CHECK(ssol_solve(scene, rng, 10, NULL, NULL), RES_BAD_ARG);
+ CHECK(ssol_solve(NULL, rng, 10, 0, NULL, &estimator), RES_BAD_ARG);
+ CHECK(ssol_solve(scene, NULL, 10, 0, NULL, &estimator), RES_BAD_ARG);
+ CHECK(ssol_solve(scene, rng, 0, 0, NULL, &estimator), RES_BAD_ARG);
+ CHECK(ssol_solve(scene, rng, 10, 0, NULL, &estimator), RES_BAD_ARG);
+ CHECK(ssol_solve(scene, rng, 10, 0, NULL, NULL), RES_BAD_ARG);
/* No geometry */
- CHECK(ssol_solve(scene, rng, 10, NULL, &estimator), RES_BAD_ARG);
+ CHECK(ssol_solve(scene, rng, 10, 0, NULL, &estimator), RES_BAD_ARG);
/* Create scene content */
CHECK(ssol_shape_create_mesh(dev, &dummy), RES_OK);
@@ -163,7 +165,40 @@ main(int argc, char** argv)
CHECK(ssol_instance_set_receiver(target, SSOL_FRONT, 0), RES_OK);
CHECK(ssol_scene_attach_instance(scene, target), RES_OK);
- CHECK(ssol_solve(scene, rng, 1, NULL, &estimator), RES_OK);
+ CHECK(ssol_solve
+ (scene, rng, 1, &SSOL_PATH_TRACKER_DEFAULT, NULL, &estimator), RES_OK);
+
+ CHECK(ssol_estimator_get_tracked_paths_count(NULL, NULL), RES_BAD_ARG);
+ CHECK(ssol_estimator_get_tracked_paths_count(estimator, NULL), RES_BAD_ARG);
+ CHECK(ssol_estimator_get_tracked_paths_count(NULL, &count), RES_BAD_ARG);
+ CHECK(ssol_estimator_get_tracked_paths_count(estimator, &count), RES_OK);
+ CHECK(count, 1);
+
+ CHECK(ssol_estimator_get_tracked_path(NULL, count, NULL), RES_BAD_ARG);
+ CHECK(ssol_estimator_get_tracked_path(estimator, count, NULL), RES_BAD_ARG);
+ CHECK(ssol_estimator_get_tracked_path(NULL, 0, NULL), RES_BAD_ARG);
+ CHECK(ssol_estimator_get_tracked_path(estimator, 0, NULL), RES_BAD_ARG);
+ CHECK(ssol_estimator_get_tracked_path(NULL, count, &path), RES_BAD_ARG);
+ CHECK(ssol_estimator_get_tracked_path(estimator, count, &path), RES_BAD_ARG);
+ CHECK(ssol_estimator_get_tracked_path(NULL, 0, &path), RES_BAD_ARG);
+ CHECK(ssol_estimator_get_tracked_path(estimator, 0, &path), RES_OK);
+
+ CHECK(ssol_path_get_vertices_count(NULL, NULL), RES_BAD_ARG);
+ CHECK(ssol_path_get_vertices_count(&path, NULL), RES_BAD_ARG);
+ CHECK(ssol_path_get_vertices_count(NULL, &count), RES_BAD_ARG);
+ CHECK(ssol_path_get_vertices_count(&path, &count), RES_OK);
+ NCHECK(count, 0);
+
+ CHECK(ssol_path_get_vertex(NULL, count, NULL), RES_BAD_ARG);
+ CHECK(ssol_path_get_vertex(&path, count, NULL), RES_BAD_ARG);
+ CHECK(ssol_path_get_vertex(NULL, 0, NULL), RES_BAD_ARG);
+ CHECK(ssol_path_get_vertex(&path, 0, NULL), RES_BAD_ARG);
+ CHECK(ssol_path_get_vertex(NULL, count, &vertex), RES_BAD_ARG);
+ CHECK(ssol_path_get_vertex(&path, count, &vertex), RES_BAD_ARG);
+ CHECK(ssol_path_get_vertex(NULL, 0, &vertex), RES_BAD_ARG);
+ FOR_EACH(i, 0, count) {
+ CHECK(ssol_path_get_vertex(&path, i, &vertex), RES_OK);
+ }
CHECK(ssol_estimator_get_sampled_area(NULL, NULL), RES_BAD_ARG);
CHECK(ssol_estimator_get_sampled_area(estimator, NULL), RES_BAD_ARG);
@@ -208,7 +243,7 @@ main(int argc, char** argv)
CHECK(ssol_instance_sample(target, 0), RES_OK);
CHECK(ssol_instance_sample(secondary, 0), RES_OK);
CHECK(ssol_instance_sample(heliostat, 0), RES_OK);
- CHECK(ssol_solve(scene, rng, 10, NULL, &estimator), RES_BAD_ARG);
+ CHECK(ssol_solve(scene, rng, 10, 0, NULL, &estimator), RES_BAD_ARG);
CHECK(ssol_estimator_get_mc_sampled(estimator, heliostat, &sampled), RES_BAD_ARG);
CHECK(ssol_instance_sample(target, 1), RES_OK);
@@ -217,7 +252,7 @@ main(int argc, char** argv)
/* No attached sun */
CHECK(ssol_scene_detach_sun(scene, sun), RES_OK);
- CHECK(ssol_solve(scene, rng, 10, NULL, &estimator), RES_BAD_ARG);
+ CHECK(ssol_solve(scene, rng, 10, 0, NULL, &estimator), RES_BAD_ARG);
CHECK(ssol_sun_ref_put(sun), RES_OK);
/* Sun with no spectrum */
@@ -225,7 +260,7 @@ main(int argc, char** argv)
CHECK(ssol_sun_set_direction(sun, d3(dir, 1, 0, -1)), RES_OK);
CHECK(ssol_sun_set_dni(sun, 1000), RES_OK);
CHECK(ssol_scene_attach_sun(scene, sun), RES_OK);
- CHECK(ssol_solve(scene, rng, 10, NULL, &estimator), RES_BAD_ARG);
+ CHECK(ssol_solve(scene, rng, 10, 0, NULL, &estimator), RES_BAD_ARG);
CHECK(ssol_scene_detach_sun(scene, sun), RES_OK);
CHECK(ssol_sun_ref_put(sun), RES_OK);
@@ -234,14 +269,14 @@ main(int argc, char** argv)
CHECK(ssol_sun_set_direction(sun, d3(dir, 1, 0, -1)), RES_OK);
CHECK(ssol_sun_set_spectrum(sun, spectrum), RES_OK);
CHECK(ssol_scene_attach_sun(scene, sun), RES_OK);
- CHECK(ssol_solve(scene, rng, 10, NULL, &estimator), RES_BAD_ARG);
+ CHECK(ssol_solve(scene, rng, 10, 0, NULL, &estimator), RES_BAD_ARG);
CHECK(ssol_sun_set_dni(sun, 1000), RES_OK);
/* No receiver in scene */
CHECK(ssol_instance_set_receiver(heliostat, 0, 0), RES_OK);
CHECK(ssol_instance_set_receiver(secondary, 0, 0), RES_OK);
CHECK(ssol_instance_set_receiver(target, 0, 0), RES_OK);
- CHECK(ssol_solve(scene, rng, 10, NULL, &estimator), RES_OK);
+ CHECK(ssol_solve(scene, rng, 10, 0, NULL, &estimator), RES_OK);
CHECK(ssol_instance_set_receiver(heliostat, SSOL_FRONT, 0), RES_OK);
CHECK(ssol_instance_set_receiver(secondary, SSOL_FRONT, 0), RES_OK);
CHECK(ssol_instance_set_receiver(target, SSOL_FRONT, 0), RES_OK);
@@ -256,7 +291,7 @@ main(int argc, char** argv)
CHECK(ssol_atmosphere_create_uniform(dev, &atm), RES_OK);
CHECK(ssol_atmosphere_set_uniform_absorption(atm, abs), RES_OK);
CHECK(ssol_scene_attach_atmosphere(scene, atm), RES_OK);
- CHECK(ssol_solve(scene, rng, 10, NULL, &estimator), RES_BAD_ARG);
+ CHECK(ssol_solve(scene, rng, 10, 0, NULL, &estimator), RES_BAD_ARG);
CHECK(ssol_scene_detach_atmosphere(scene, atm), RES_OK);
CHECK(ssol_spectrum_ref_put(abs), RES_OK);
CHECK(ssol_atmosphere_ref_put(atm), RES_OK);
@@ -266,7 +301,7 @@ main(int argc, char** argv)
#define N__ 10000
#define GET_MC_RCV ssol_estimator_get_mc_receiver
#define GET_MC_GLOBAL ssol_estimator_get_mc_global
- CHECK(ssol_solve(scene, rng, N__, tmp, &estimator), RES_OK);
+ CHECK(ssol_solve(scene, rng, N__, 0, tmp, &estimator), RES_OK);
CHECK(ssol_instance_get_id(target, &r_id), RES_OK);
CHECK(ssol_estimator_get_realisation_count(estimator, &count), RES_OK);
CHECK(count, N__);
@@ -315,7 +350,7 @@ main(int argc, char** argv)
CHECK(ssol_instance_sample(secondary, 0), RES_OK);
NCHECK(tmp = tmpfile(), 0);
- CHECK(ssol_solve(scene, rng, N__, tmp, &estimator), RES_OK);
+ CHECK(ssol_solve(scene, rng, N__, 0, tmp, &estimator), RES_OK);
CHECK(ssol_estimator_get_realisation_count(estimator, &count), RES_OK);
CHECK(count, N__);
CHECK(pp_sum(tmp, (int32_t)r_id, count, &m, &std), RES_OK);
@@ -348,7 +383,7 @@ main(int argc, char** argv)
CHECK(ssol_scene_attach_atmosphere(scene, atm), RES_OK);
NCHECK(tmp = tmpfile(), 0);
- CHECK(ssol_solve(scene, rng, N__, tmp, &estimator), RES_OK);
+ CHECK(ssol_solve(scene, rng, N__, 0, tmp, &estimator), RES_OK);
CHECK(ssol_estimator_get_realisation_count(estimator, &count), RES_OK);
CHECK(count, N__);
CHECK(pp_sum(tmp, (int32_t)r_id, count, &m, &std), RES_OK);
@@ -399,7 +434,7 @@ main(int argc, char** argv)
CHECK(ssol_instance_set_receiver(target, SSOL_FRONT, 1), RES_OK);
NCHECK(tmp = tmpfile(), 0);
- CHECK(ssol_solve(scene, rng, N__, tmp, &estimator), RES_OK);
+ CHECK(ssol_solve(scene, rng, N__, 0, tmp, &estimator), RES_OK);
CHECK(ssol_estimator_get_realisation_count(estimator, &count), RES_OK);
CHECK(count, N__);
CHECK(pp_sum(tmp, (int32_t)r_id, count, &a_m, &a_std), RES_OK);
@@ -499,7 +534,7 @@ main(int argc, char** argv)
desc.count = 2;
CHECK(ssol_spectrum_setup(abs, get_wlen, 2, &desc), RES_OK);
NCHECK(tmp = tmpfile(), 0);
- CHECK(ssol_solve(scene, rng, N__, tmp, &estimator), RES_OK);
+ CHECK(ssol_solve(scene, rng, N__, 0, tmp, &estimator), RES_OK);
CHECK(ssol_estimator_get_realisation_count(estimator, &count), RES_OK);
CHECK(count, N__);
CHECK(pp_sum(tmp, (int32_t)r_id, count, &m, &std), RES_OK);
diff --git a/src/test_ssol_solver2.c b/src/test_ssol_solver2.c
@@ -180,7 +180,7 @@ main(int argc, char** argv)
NCHECK(tmp = tmpfile(), 0);
#define N__ 10000
#define GET_MC_RCV ssol_estimator_get_mc_receiver
- CHECK(ssol_solve(scene, rng, N__, tmp, &estimator), RES_OK);
+ CHECK(ssol_solve(scene, rng, N__, 0, tmp, &estimator), RES_OK);
CHECK(ssol_instance_get_id(target, &r_id), RES_OK);
CHECK(ssol_estimator_get_realisation_count(estimator, &count), RES_OK);
CHECK(count, N__);
diff --git a/src/test_ssol_solver2b.c b/src/test_ssol_solver2b.c
@@ -183,7 +183,7 @@ main(int argc, char** argv)
NCHECK(tmp = tmpfile(), 0);
#define N__ 50000
- CHECK(ssol_solve(scene, rng, N__, tmp, &estimator), RES_OK);
+ CHECK(ssol_solve(scene, rng, N__, 0, tmp, &estimator), RES_OK);
CHECK(ssol_instance_get_id(target, &r_id), RES_OK);
CHECK(ssol_estimator_get_realisation_count(estimator, &count), RES_OK);
CHECK(count, N__);
diff --git a/src/test_ssol_solver3.c b/src/test_ssol_solver3.c
@@ -137,7 +137,7 @@ main(int argc, char** argv)
NCHECK(tmp = tmpfile(), 0);
#define N__ 20000
- CHECK(ssol_solve(scene, rng, N__, tmp, &estimator), RES_OK);
+ CHECK(ssol_solve(scene, rng, N__, 0, tmp, &estimator), RES_OK);
CHECK(ssol_instance_get_id(target, &r_id), RES_OK);
CHECK(ssol_estimator_get_realisation_count(estimator, &count), RES_OK);
CHECK(count, N__);
diff --git a/src/test_ssol_solver4.c b/src/test_ssol_solver4.c
@@ -146,7 +146,7 @@ main(int argc, char** argv)
NCHECK(tmp = tmpfile(), 0);
#define N__ 10000
#define GET_MC_RCV ssol_estimator_get_mc_receiver
- CHECK(ssol_solve(scene, rng, N__, tmp, &estimator), RES_OK);
+ CHECK(ssol_solve(scene, rng, N__, 0, tmp, &estimator), RES_OK);
CHECK(ssol_instance_get_id(target1, &r_id1), RES_OK);
CHECK(ssol_instance_get_id(target2, &r_id2), RES_OK);
CHECK(ssol_estimator_get_realisation_count(estimator, &count), RES_OK);
diff --git a/src/test_ssol_solver5.c b/src/test_ssol_solver5.c
@@ -138,7 +138,7 @@ main(int argc, char** argv)
NCHECK(tmp = tmpfile(), 0);
#define N__ 10000
- CHECK(ssol_solve(scene, rng, N__, tmp, &estimator), RES_OK);
+ CHECK(ssol_solve(scene, rng, N__, 0, tmp, &estimator), RES_OK);
CHECK(ssol_instance_get_id(target, &r_id), RES_OK);
CHECK(ssol_estimator_get_realisation_count(estimator, &count), RES_OK);
CHECK(count, N__);
diff --git a/src/test_ssol_solver6.c b/src/test_ssol_solver6.c
@@ -178,7 +178,7 @@ main(int argc, char** argv)
NCHECK(tmp = tmpfile(), 0);
#define N__ 10000
#define GET_MC_RCV ssol_estimator_get_mc_receiver
- CHECK(ssol_solve(scene, rng, N__, tmp, &estimator), RES_OK);
+ CHECK(ssol_solve(scene, rng, N__, 0, tmp, &estimator), RES_OK);
CHECK(fclose(tmp), 0);
CHECK(ssol_estimator_get_mc_global(estimator, &mc_global), RES_OK);
printf("Shadows = %g +/- %g\n", mc_global.shadowed.E, mc_global.shadowed.SE);
diff --git a/src/test_ssol_solver7.c b/src/test_ssol_solver7.c
@@ -198,7 +198,7 @@ main(int argc, char** argv)
#define GET_MC_RCV ssol_estimator_get_mc_receiver
NCHECK(tmp = tmpfile(), 0);
- CHECK(ssol_solve(scene, rng, N__, tmp, &estimator), RES_OK);
+ CHECK(ssol_solve(scene, rng, N__, 0, tmp, &estimator), RES_OK);
CHECK(fclose(tmp), 0);
printf("Total = %g\n", TOTAL);
diff --git a/src/test_ssol_solver8.c b/src/test_ssol_solver8.c
@@ -0,0 +1,186 @@
+/* Copyright (C) CNRS 2016-2017
+*
+* This program is free software: you can redistribute it and/or modify
+* it under the terms of the GNU General Public License as published by
+* the Free Software Foundation, either version 3 of the License, or
+* (at your option) any later version.
+*
+* This program is distributed in the hope that it will be useful,
+* but WITHOUT ANY WARRANTY; without even the implied warranty of
+* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+* GNU General Public License for more details.
+*
+* You should have received a copy of the GNU General Public License
+* along with this program. If not, see <http://www.gnu.org/licenses/>. */
+
+#include "ssol.h"
+#include "test_ssol_utils.h"
+#include "test_ssol_materials.h"
+
+#define PLANE_NAME TARGET
+#define HALF_X 10
+#define HALF_Y 10
+#include "test_ssol_rect_geometry.h"
+
+#define X_SZ 10
+#define Y_SZ 4
+#define POLYGON_NAME POLY
+#define HALF_X (X_SZ / 2)
+STATIC_ASSERT((HALF_X * 2 == X_SZ), ONLY_ENVEN_VALUES_FOR_X_SZ);
+#define Y_MIN 0
+#define Y_MAX Y_SZ
+#include "test_ssol_rect2D_geometry.h"
+
+#include <rsys/double33.h>
+
+#include <star/s3d.h>
+#include <star/ssp.h>
+
+static void
+get_wlen(const size_t i, double* wlen, double* data, void* ctx)
+{
+ double wavelengths[3] = { 1, 2, 3 };
+ double intensities[3] = { 1, 0.8, 1 };
+ CHECK(i < 3, 1);
+ (void) ctx;
+ *wlen = wavelengths[i];
+ *data = intensities[i];
+}
+
+int
+main(int argc, char** argv)
+{
+ struct mem_allocator allocator;
+ struct ssol_device* dev;
+ struct ssp_rng* rng;
+ struct ssol_scene* scene;
+ struct ssol_shape* square;
+ struct ssol_vertex_data attribs[1] = { SSOL_VERTEX_DATA_NULL__ };
+ struct ssol_shape* quad_square;
+ struct ssol_carving carving = SSOL_CARVING_NULL;
+ struct ssol_quadric quadric = SSOL_QUADRIC_DEFAULT;
+ struct ssol_punched_surface punched = SSOL_PUNCHED_SURFACE_NULL;
+ struct ssol_material* m_mtl;
+ struct ssol_material* v_mtl;
+ struct ssol_mirror_shader shader = SSOL_MIRROR_SHADER_NULL;
+ struct ssol_object* m_object;
+ struct ssol_object* t_object;
+ struct ssol_instance* heliostat;
+ struct ssol_instance* target;
+ struct ssol_sun* sun;
+ struct ssol_spectrum* spectrum;
+ struct ssol_estimator* estimator;
+ struct ssol_mc_global mc_global;
+ struct ssol_mc_receiver mc_rcv;
+ double dir[3];
+ double transform[12]; /* 3x4 column major matrix */
+ size_t count;
+ FILE* tmp;
+ uint32_t r_id;
+
+ (void) argc, (void) argv;
+ d3_splat(transform + 9, 0);
+ d33_rotation_pitch(transform, PI); /* flip faces: invert normal */
+ transform[11] = 4; /* set it just above the parabolic cylinder */
+
+ mem_init_proxy_allocator(&allocator, &mem_default_allocator);
+
+ CHECK(ssol_device_create
+ (NULL, &allocator, SSOL_NTHREADS_DEFAULT, 0, &dev), RES_OK);
+
+#define DNI 1000
+ CHECK(ssp_rng_create(&allocator, &ssp_rng_threefry, &rng), RES_OK);
+ CHECK(ssol_spectrum_create(dev, &spectrum), RES_OK);
+ CHECK(ssol_spectrum_setup(spectrum, get_wlen, 3, NULL), RES_OK);
+ CHECK(ssol_sun_create_directional(dev, &sun), RES_OK);
+ CHECK(ssol_sun_set_direction(sun, d3(dir, 0, 1, -1)), RES_OK);
+ CHECK(ssol_sun_set_spectrum(sun, spectrum), RES_OK);
+ CHECK(ssol_sun_set_dni(sun, DNI), RES_OK);
+ CHECK(ssol_scene_create(dev, &scene), RES_OK);
+ CHECK(ssol_scene_attach_sun(scene, sun), RES_OK);
+
+ /* Create scene content */
+
+ CHECK(ssol_shape_create_mesh(dev, &square), RES_OK);
+ attribs[0].usage = SSOL_POSITION;
+ attribs[0].get = get_position;
+ CHECK(ssol_mesh_setup(square, TARGET_NTRIS__, get_ids,
+ TARGET_NVERTS__, attribs, 1, (void*) &TARGET_DESC__), RES_OK);
+
+ CHECK(ssol_shape_create_punched_surface(dev, &quad_square), RES_OK);
+ carving.get = get_polygon_vertices;
+ carving.operation = SSOL_AND;
+ carving.nb_vertices = POLY_NVERTS__;
+ carving.context = &POLY_EDGES__;
+ quadric.type = SSOL_QUADRIC_PARABOLIC_CYLINDER;
+ quadric.data.parabol.focal = 1;
+ punched.nb_carvings = 1;
+ punched.quadric = &quadric;
+ punched.carvings = &carving;
+ CHECK(ssol_punched_surface_setup(quad_square, &punched), RES_OK);
+
+ CHECK(ssol_material_create_mirror(dev, &m_mtl), RES_OK);
+ shader.normal = get_shader_normal;
+ shader.reflectivity = get_shader_reflectivity;
+ shader.roughness = get_shader_roughness;
+ CHECK(ssol_mirror_set_shader(m_mtl, &shader), RES_OK);
+ CHECK(ssol_material_create_virtual(dev, &v_mtl), RES_OK);
+
+ CHECK(ssol_object_create(dev, &m_object), RES_OK);
+ CHECK(ssol_object_add_shaded_shape(m_object, quad_square, m_mtl, m_mtl), RES_OK);
+ CHECK(ssol_object_instantiate(m_object, &heliostat), RES_OK);
+ CHECK(ssol_scene_attach_instance(scene, heliostat), RES_OK);
+
+ CHECK(ssol_object_create(dev, &t_object), RES_OK);
+ CHECK(ssol_object_add_shaded_shape(t_object, square, v_mtl, v_mtl), RES_OK);
+ CHECK(ssol_object_instantiate(t_object, &target), RES_OK);
+ CHECK(ssol_instance_set_transform(target, transform), RES_OK);
+ CHECK(ssol_instance_set_receiver(target, SSOL_FRONT, 0), RES_OK);
+ CHECK(ssol_instance_sample(target, 0), RES_OK);
+ CHECK(ssol_scene_attach_instance(scene, target), RES_OK);
+
+ NCHECK(tmp = tmpfile(), 0);
+#define N__ 100000
+#define GET_MC_RCV ssol_estimator_get_mc_receiver
+ CHECK(ssol_solve(scene, rng, N__, 0, tmp, &estimator), RES_OK);
+ CHECK(ssol_instance_get_id(target, &r_id), RES_OK);
+ CHECK(ssol_estimator_get_count(estimator, &count), RES_OK);
+ CHECK(count, N__);
+ CHECK(fclose(tmp), 0);
+ CHECK(ssol_estimator_get_failed_count(estimator, &count), RES_OK);
+ CHECK(count, 0);
+#define S (sqrt(2) * X_SZ * Y_SZ)
+ CHECK(ssol_estimator_get_mc_global(estimator, &mc_global), RES_OK);
+ printf("Cos = %g +/- %g\n", mc_global.cos_loss.E, mc_global.cos_loss.SE);
+ printf("Shadows = %g +/- %g\n", mc_global.shadowed.E, mc_global.shadowed.SE);
+ printf("Missing = %g +/- %g\n", mc_global.missing.E, mc_global.missing.SE);
+ CHECK(eq_eps(mc_global.cos_loss.E, 0, 1e-4), 1);
+ CHECK(eq_eps(mc_global.shadowed.E, 0, 1e-4), 1);
+ CHECK(eq_eps(mc_global.missing.E, 0, 1e-4), 1);
+ CHECK(GET_MC_RCV(estimator, target, SSOL_FRONT, &mc_rcv), RES_OK);
+ printf("Ir(target1) = %g +/- %g\n",
+ mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE);
+ CHECK(eq_eps(mc_rcv.integrated_irradiance.E, S * DNI, 2 * mc_rcv.integrated_irradiance.SE), 1);
+
+ /* Free data */
+ CHECK(ssol_instance_ref_put(heliostat), RES_OK);
+ CHECK(ssol_instance_ref_put(target), RES_OK);
+ CHECK(ssol_object_ref_put(m_object), RES_OK);
+ CHECK(ssol_object_ref_put(t_object), RES_OK);
+ CHECK(ssol_shape_ref_put(square), RES_OK);
+ CHECK(ssol_shape_ref_put(quad_square), RES_OK);
+ CHECK(ssol_material_ref_put(m_mtl), RES_OK);
+ CHECK(ssol_material_ref_put(v_mtl), RES_OK);
+ CHECK(ssol_estimator_ref_put(estimator), RES_OK);
+ CHECK(ssol_device_ref_put(dev), RES_OK);
+ CHECK(ssol_scene_ref_put(scene), RES_OK);
+ CHECK(ssp_rng_ref_put(rng), RES_OK);
+ CHECK(ssol_spectrum_ref_put(spectrum), RES_OK);
+ CHECK(ssol_sun_ref_put(sun), RES_OK);
+
+ check_memory_allocator(&allocator);
+ mem_shutdown_proxy_allocator(&allocator);
+ CHECK(mem_allocated_size(), 0);
+
+ return 0;
+}
diff --git a/src/test_ssol_solver9.c b/src/test_ssol_solver9.c
@@ -0,0 +1,184 @@
+/* Copyright (C) CNRS 2016-2017
+*
+* This program is free software: you can redistribute it and/or modify
+* it under the terms of the GNU General Public License as published by
+* the Free Software Foundation, either version 3 of the License, or
+* (at your option) any later version.
+*
+* This program is distributed in the hope that it will be useful,
+* but WITHOUT ANY WARRANTY; without even the implied warranty of
+* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+* GNU General Public License for more details.
+*
+* You should have received a copy of the GNU General Public License
+* along with this program. If not, see <http://www.gnu.org/licenses/>. */
+
+#include "ssol.h"
+#include "test_ssol_utils.h"
+#include "test_ssol_materials.h"
+
+#include <rsys/math.h>
+
+#define TGT_X 6
+#define TGT_Y 10
+#define PLANE_NAME TARGET
+#define HALF_X (TGT_X / 2)
+#define HALF_Y (TGT_Y / 2)
+STATIC_ASSERT((HALF_X * 2 == TGT_X), ONLY_ENVEN_VALUES_FOR_TGT_X);
+STATIC_ASSERT((HALF_Y * 2 == TGT_Y), ONLY_ENVEN_VALUES_FOR_TGT_Y);
+#include "test_ssol_rect_geometry.h"
+
+#define SZ MMAX(TGT_X, TGT_Y)
+#define CUBE_NAME CUBE
+#define HALF_X (SZ / 2)
+#define HALF_Y (SZ / 2)
+#define HALF_Z (SZ / 2)
+STATIC_ASSERT((HALF_X * 2 == SZ), ONLY_ENVEN_VALUES_FOR_SZ);
+#include "test_ssol_cube_geometry.h"
+
+#include <rsys/double33.h>
+
+#include <star/s3d.h>
+#include <star/ssp.h>
+
+static void
+get_wlen(const size_t i, double* wlen, double* data, void* ctx)
+{
+ double wavelengths[3] = { 1, 2, 3 };
+ double intensities[3] = { 1, 0.8, 1 };
+ CHECK(i < 3, 1);
+ (void) ctx;
+ *wlen = wavelengths[i];
+ *data = intensities[i];
+}
+
+int
+main(int argc, char** argv)
+{
+ struct mem_allocator allocator;
+ struct ssol_device* dev;
+ struct ssp_rng* rng;
+ struct ssol_scene* scene;
+ struct ssol_shape* square;
+ struct ssol_vertex_data attribs[1] = { SSOL_VERTEX_DATA_NULL__ };
+ struct ssol_shape* cube;
+ struct ssol_material* m_mtl;
+ struct ssol_material* v_mtl;
+ struct ssol_mirror_shader shader = SSOL_MIRROR_SHADER_NULL;
+ struct ssol_object* m_object;
+ struct ssol_object* t_object;
+ struct ssol_instance* heliostat;
+ struct ssol_instance* target;
+ struct ssol_sun* sun;
+ struct ssol_spectrum* spectrum;
+ struct ssol_estimator* estimator;
+ struct ssol_mc_global mc_global;
+ struct ssol_mc_receiver mc_rcv;
+ double dir[3];
+ double transform[12]; /* 3x4 column major matrix */
+ size_t count;
+ FILE* tmp;
+ uint32_t r_id;
+
+ (void) argc, (void) argv;
+ d3_splat(transform + 9, 0);
+ d33_rotation_pitch(transform, PI); /* flip faces: invert normal */
+ transform[11] = 6; /* set it above the cube */
+
+ mem_init_proxy_allocator(&allocator, &mem_default_allocator);
+
+ CHECK(ssol_device_create
+ (NULL, &allocator, SSOL_NTHREADS_DEFAULT, 0, &dev), RES_OK);
+
+#define DNI 1000
+ CHECK(ssp_rng_create(&allocator, &ssp_rng_threefry, &rng), RES_OK);
+ CHECK(ssol_spectrum_create(dev, &spectrum), RES_OK);
+ CHECK(ssol_spectrum_setup(spectrum, get_wlen, 3, NULL), RES_OK);
+ CHECK(ssol_sun_create_directional(dev, &sun), RES_OK);
+ CHECK(ssol_sun_set_direction(sun, d3(dir, 0, 0, -1)), RES_OK);
+ CHECK(ssol_sun_set_spectrum(sun, spectrum), RES_OK);
+ CHECK(ssol_sun_set_dni(sun, DNI), RES_OK);
+ CHECK(ssol_scene_create(dev, &scene), RES_OK);
+ CHECK(ssol_scene_attach_sun(scene, sun), RES_OK);
+
+ /* Create scene content */
+
+ CHECK(ssol_shape_create_mesh(dev, &square), RES_OK);
+ attribs[0].usage = SSOL_POSITION;
+ attribs[0].get = get_position;
+ CHECK(ssol_mesh_setup(square, TARGET_NTRIS__, get_ids,
+ TARGET_NVERTS__, attribs, 1, (void*) &TARGET_DESC__), RES_OK);
+
+ CHECK(ssol_shape_create_mesh(dev, &cube), RES_OK);
+ CHECK(ssol_mesh_setup(cube, CUBE_NTRIS__, get_ids,
+ CUBE_NVERTS__, attribs, 1, (void*) &CUBE_DESC__), RES_OK);
+
+ CHECK(ssol_material_create_mirror(dev, &m_mtl), RES_OK);
+ shader.normal = get_shader_normal;
+ shader.reflectivity = get_shader_reflectivity;
+ shader.roughness = get_shader_roughness;
+ CHECK(ssol_mirror_set_shader(m_mtl, &shader), RES_OK);
+ CHECK(ssol_material_create_virtual(dev, &v_mtl), RES_OK);
+
+ CHECK(ssol_object_create(dev, &m_object), RES_OK);
+ CHECK(ssol_object_add_shaded_shape(m_object, cube, m_mtl, m_mtl), RES_OK);
+ CHECK(ssol_object_instantiate(m_object, &heliostat), RES_OK);
+ CHECK(ssol_scene_attach_instance(scene, heliostat), RES_OK);
+
+ CHECK(ssol_object_create(dev, &t_object), RES_OK);
+ CHECK(ssol_object_add_shaded_shape(t_object, square, v_mtl, v_mtl), RES_OK);
+ CHECK(ssol_object_instantiate(t_object, &target), RES_OK);
+ CHECK(ssol_instance_set_transform(target, transform), RES_OK);
+ CHECK(ssol_instance_set_receiver(target, SSOL_FRONT, 0), RES_OK);
+ CHECK(ssol_instance_sample(target, 0), RES_OK);
+ CHECK(ssol_scene_attach_instance(scene, target), RES_OK);
+
+ NCHECK(tmp = tmpfile(), 0);
+#define N__ 100000
+#define GET_MC_RCV ssol_estimator_get_mc_receiver
+ CHECK(ssol_solve(scene, rng, N__, 0, tmp, &estimator), RES_OK);
+ CHECK(ssol_instance_get_id(target, &r_id), RES_OK);
+ CHECK(ssol_estimator_get_count(estimator, &count), RES_OK);
+ CHECK(count, N__);
+ CHECK(fclose(tmp), 0);
+ CHECK(ssol_estimator_get_failed_count(estimator, &count), RES_OK);
+ CHECK(count, 0);
+#define DNI_TGT_S (DNI * TGT_X * TGT_Y)
+#define DNI_S (DNI * SZ * SZ)
+ CHECK(ssol_estimator_get_mc_global(estimator, &mc_global), RES_OK);
+ printf("Cos = %g +/- %g\n", mc_global.cos_loss.E, mc_global.cos_loss.SE);
+ printf("Shadows = %g +/- %g\n", mc_global.shadowed.E, mc_global.shadowed.SE);
+ printf("Missing = %g +/- %g\n", mc_global.missing.E, mc_global.missing.SE);
+ /* CHECK(eq_eps(mc_global.cos_loss.E, 4 * DNI_S, DNI_S * 1e-4), 1); */
+ CHECK(eq_eps(mc_global.shadowed.E, DNI_S, 3 * mc_global.shadowed.SE), 1);
+ CHECK(eq_eps(mc_global.missing.E,
+ MMAX(DNI_S, DNI_TGT_S) - MMIN(DNI_S, DNI_TGT_S),
+ 3 * mc_global.missing.SE), 1);
+ CHECK(GET_MC_RCV(estimator, target, SSOL_FRONT, &mc_rcv), RES_OK);
+ printf("Ir(target1) = %g +/- %g\n",
+ mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE);
+ CHECK(eq_eps(mc_rcv.integrated_irradiance.E, MMIN(DNI_S, DNI_TGT_S),
+ 3 * mc_rcv.integrated_irradiance.SE), 1);
+
+ /* Free data */
+ CHECK(ssol_instance_ref_put(heliostat), RES_OK);
+ CHECK(ssol_instance_ref_put(target), RES_OK);
+ CHECK(ssol_object_ref_put(m_object), RES_OK);
+ CHECK(ssol_object_ref_put(t_object), RES_OK);
+ CHECK(ssol_shape_ref_put(square), RES_OK);
+ CHECK(ssol_shape_ref_put(cube), RES_OK);
+ CHECK(ssol_material_ref_put(m_mtl), RES_OK);
+ CHECK(ssol_material_ref_put(v_mtl), RES_OK);
+ CHECK(ssol_estimator_ref_put(estimator), RES_OK);
+ CHECK(ssol_device_ref_put(dev), RES_OK);
+ CHECK(ssol_scene_ref_put(scene), RES_OK);
+ CHECK(ssp_rng_ref_put(rng), RES_OK);
+ CHECK(ssol_spectrum_ref_put(spectrum), RES_OK);
+ CHECK(ssol_sun_ref_put(sun), RES_OK);
+
+ check_memory_allocator(&allocator);
+ mem_shutdown_proxy_allocator(&allocator);
+ CHECK(mem_allocated_size(), 0);
+
+ return 0;
+}
