solstice-solver

Solver library of the solstice app
git clone git://git.meso-star.com/solstice-solver.git
Log | Files | Refs | README | LICENSE

commit 45d77042123d07e5dc6f18d88c69734118da96ae
parent 7bf7d7f4bb3dae6effaf2d94b56f237bbf3724cc
Author: Christophe Coustet <christophe.coustet@meso-star.com>
Date:   Fri, 31 Mar 2017 11:06:51 +0200

Merge remote-tracking branch 'origin/feature_outputs' into develop

Diffstat:
Msrc/ssol.h | 31+++++++++++++++++++++++++++----
Msrc/ssol_estimator.c | 45++++++++++++++++++++++++++++++++++++++++-----
Msrc/ssol_estimator_c.h | 38++++++++++++++++----------------------
Msrc/ssol_instance.c | 2+-
Msrc/ssol_mc_receiver.c | 6+++---
Msrc/ssol_object.c | 9++++++---
Msrc/ssol_object_c.h | 1+
Msrc/ssol_shape.c | 49++++++++++++++++++++++++++++++-------------------
Msrc/ssol_shape_c.h | 11++++++-----
Msrc/ssol_solver.c | 170++++++++++++++++++++++++++++++++++++++++++++++++-------------------------------
Msrc/test_ssol_instance.c | 6+++++-
Msrc/test_ssol_object.c | 5+++++
Msrc/test_ssol_solver1.c | 73++++++++++++++++++++++++++++++++++++++++++++++++-------------------------
Msrc/test_ssol_solver2.c | 11+++++++----
Msrc/test_ssol_solver2b.c | 9++++++---
Msrc/test_ssol_solver3.c | 12+++++++++---
Msrc/test_ssol_solver4.c | 19+++++++++++--------
Msrc/test_ssol_solver5.c | 8+++++---
Msrc/test_ssol_solver6.c | 4++++
Msrc/test_ssol_solver7.c | 29+++++++++++++++--------------
Msrc/test_ssol_solver8.c | 8++++----
Msrc/test_ssol_solver9.c | 8+++-----
22 files changed, 355 insertions(+), 199 deletions(-)

diff --git a/src/ssol.h b/src/ssol.h @@ -347,22 +347,28 @@ struct ssol_mc_result { static const struct ssol_mc_result SSOL_MC_RESULT_NULL = SSOL_MC_RESULT_NULL__; struct ssol_mc_global { - struct ssol_mc_result cos_loss; /* In W */ + struct ssol_mc_result cos_factor; /* [0 1] */ + struct ssol_mc_result absorbed; /* In W */ struct ssol_mc_result shadowed; /* In W */ struct ssol_mc_result missing; /* In W */ + struct ssol_mc_result atmosphere; /* In W */ + struct ssol_mc_result reflectivity; /* In W */ }; #define SSOL_MC_GLOBAL_NULL__ { \ SSOL_MC_RESULT_NULL__, \ SSOL_MC_RESULT_NULL__, \ + SSOL_MC_RESULT_NULL__, \ + SSOL_MC_RESULT_NULL__, \ + SSOL_MC_RESULT_NULL__, \ SSOL_MC_RESULT_NULL__ \ } static const struct ssol_mc_global SSOL_MC_GLOBAL_NULL = SSOL_MC_GLOBAL_NULL__; struct ssol_mc_receiver { struct ssol_mc_result integrated_irradiance; /* In W */ + struct ssol_mc_result integrated_absorbed_irradiance; /* In W */ struct ssol_mc_result absorptivity_loss; /* In W */ struct ssol_mc_result reflectivity_loss; /* In W */ - struct ssol_mc_result cos_loss; /* In W TODO remove this */ /* Internal data */ size_t N__; @@ -388,11 +394,17 @@ struct ssol_mc_shape { #define SSOL_MC_SHAPE_NULL__ { 0, NULL, NULL } static const struct ssol_mc_shape SSOL_MC_SHAPE_NULL = SSOL_MC_SHAPE_NULL__; +struct ssol_mc_sampled { + struct ssol_mc_result cos_factor; /* [0 1] */ + struct ssol_mc_result shadowed; + size_t nb_samples; +}; + struct ssol_mc_primitive { struct ssol_mc_result integrated_irradiance; /* In W */ + struct ssol_mc_result integrated_absorbed_irradiance; /* In W */ struct ssol_mc_result absorptivity_loss; /* In W */ struct ssol_mc_result reflectivity_loss; /* In W */ - struct ssol_mc_result cos_loss; /* In W TODO remove this */ }; #define SSOL_MC_PRIMITIVE_NULL__ { \ SSOL_MC_RESULT_NULL__, \ @@ -1014,7 +1026,7 @@ ssol_estimator_get_mc_sampled_x_receiver struct ssol_mc_receiver* rcv); SSOL_API res_T -ssol_estimator_get_count +ssol_estimator_get_realisation_count (const struct ssol_estimator* estimator, size_t* count); @@ -1029,6 +1041,17 @@ ssol_estimator_get_sampled_area (const struct ssol_estimator* estimator, double* area); +SSOL_API res_T +ssol_estimator_get_sampled_count + (const struct ssol_estimator* estimator, + size_t* count); + +SSOL_API res_T +ssol_estimator_get_mc_sampled + (struct ssol_estimator* estimator, + const struct ssol_instance* samp_instance, + struct ssol_mc_sampled* sampled); + /******************************************************************************* * Tracked paths ******************************************************************************/ diff --git a/src/ssol_estimator.c b/src/ssol_estimator.c @@ -118,9 +118,12 @@ ssol_estimator_get_mc_global global->Name.V = data->sqr_weight / N - global->Name.E*global->Name.E; \ global->Name.SE = global->Name.V > 0 ? sqrt(global->Name.V / N) : 0; \ } (void)0 - SETUP_MC_RESULT(cos_loss); + SETUP_MC_RESULT(cos_factor); + SETUP_MC_RESULT(absorbed); SETUP_MC_RESULT(shadowed); SETUP_MC_RESULT(missing); + SETUP_MC_RESULT(atmosphere); + SETUP_MC_RESULT(reflectivity); #undef SETUP_MC_RESULT return RES_OK; } @@ -145,7 +148,6 @@ ssol_estimator_get_mc_sampled_x_receiver memset(rcv, 0, sizeof(rcv[0])); - mc_samp = htable_sampled_find(&estimator->mc_sampled, &samp_instance); if(!mc_samp || !mc_samp->nb_samples) { /* The sampled instance has no MC estimation */ @@ -161,16 +163,16 @@ ssol_estimator_get_mc_sampled_x_receiver mc_rcv1 = side == SSOL_FRONT ? &mc_rcv->front : &mc_rcv->back; #define SETUP_MC_RESULT(Name) { \ - const double N = (double)mc_samp->nb_samples; \ + const double N = (double)estimator->realisation_count; \ const struct mc_data* data = &mc_rcv1->Name; \ rcv->Name.E = data->weight / N; \ rcv->Name.V = data->sqr_weight / N - rcv->Name.E*rcv->Name.E; \ rcv->Name.SE = rcv->Name.V > 0 ? sqrt(rcv->Name.V / N) : 0; \ } (void)0 SETUP_MC_RESULT(integrated_irradiance); + SETUP_MC_RESULT(integrated_absorbed_irradiance); SETUP_MC_RESULT(absorptivity_loss); SETUP_MC_RESULT(reflectivity_loss); - SETUP_MC_RESULT(cos_loss); #undef SETUP_MC_RESULT rcv->mc__ = mc_rcv1; rcv->N__ = mc_samp->nb_samples; @@ -178,7 +180,7 @@ ssol_estimator_get_mc_sampled_x_receiver } res_T -ssol_estimator_get_count +ssol_estimator_get_realisation_count (const struct ssol_estimator* estimator, size_t* count) { if (!estimator || !count) return RES_BAD_ARG; @@ -205,6 +207,39 @@ ssol_estimator_get_sampled_area } res_T +ssol_estimator_get_sampled_count + (const struct ssol_estimator* estimator, size_t* count) +{ + if (!estimator || !count) return RES_BAD_ARG; + *count = htable_sampled_size_get(&estimator->mc_sampled); + return RES_OK; +} + +res_T +ssol_estimator_get_mc_sampled + (struct ssol_estimator* estimator, + const struct ssol_instance* samp_instance, + struct ssol_mc_sampled* sampled) +{ + struct mc_sampled* mc = NULL; + if (!estimator || !samp_instance || !sampled) return RES_BAD_ARG; + mc = htable_sampled_find(&estimator->mc_sampled, &samp_instance); + if(!mc) return RES_BAD_ARG; + sampled->nb_samples = mc->nb_samples; + #define SETUP_MC_RESULT(Name) { \ + const double N = (double)estimator->realisation_count; \ + const struct mc_data* data = &mc->Name; \ + sampled->Name.E = data->weight / N; \ + sampled->Name.V = data->sqr_weight/N - sampled->Name.E*sampled->Name.E; \ + sampled->Name.SE = sampled->Name.V > 0 ? sqrt(sampled->Name.V / N) : 0; \ + } (void)0 + SETUP_MC_RESULT(cos_factor); + SETUP_MC_RESULT(shadowed); + #undef SETUP_MC_RESULT + return RES_OK; +} + +res_T ssol_estimator_get_tracked_paths_count (const struct ssol_estimator* estimator, size_t* npaths) { diff --git a/src/ssol_estimator_c.h b/src/ssol_estimator_c.h @@ -37,9 +37,9 @@ static const struct mc_data MC_DATA_NULL = MC_DATA_NULL__; #define MC_RECEIVER_DATA \ struct mc_data integrated_irradiance; /* In W */ \ + struct mc_data integrated_absorbed_irradiance; /* In W */ \ struct mc_data absorptivity_loss; /* In W */ \ - struct mc_data reflectivity_loss; /* In W */ \ - struct mc_data cos_loss; /* In W */ + struct mc_data reflectivity_loss; /* In W */ #define MC_RECEIVER_DATA_NULL__ \ MC_DATA_NULL__, \ @@ -148,9 +148,9 @@ mc_receiver_1side_init { ASSERT(mc); mc->integrated_irradiance = MC_DATA_NULL; + mc->integrated_absorbed_irradiance = MC_DATA_NULL; mc->absorptivity_loss = MC_DATA_NULL; mc->reflectivity_loss = MC_DATA_NULL; - mc->cos_loss = MC_DATA_NULL; htable_shape2mc_init(allocator, &mc->shape2mc); } @@ -167,9 +167,9 @@ mc_receiver_1side_copy { ASSERT(dst && src); dst->integrated_irradiance = src->integrated_irradiance; + dst->integrated_absorbed_irradiance = src->integrated_absorbed_irradiance; dst->absorptivity_loss = src->absorptivity_loss; dst->reflectivity_loss = src->reflectivity_loss; - dst->cos_loss = src->cos_loss; return htable_shape2mc_copy(&dst->shape2mc, &src->shape2mc); } @@ -179,9 +179,9 @@ mc_receiver_1side_copy_and_release { ASSERT(dst && src); dst->integrated_irradiance = src->integrated_irradiance; + dst->integrated_absorbed_irradiance = src->integrated_absorbed_irradiance; dst->absorptivity_loss = src->absorptivity_loss; dst->reflectivity_loss = src->reflectivity_loss; - dst->cos_loss = src->cos_loss; return htable_shape2mc_copy_and_release(&dst->shape2mc, &src->shape2mc); } @@ -278,10 +278,8 @@ mc_receiver_copy_and_release ******************************************************************************/ struct mc_sampled { /* Global data for this entity */ - struct mc_data cos_loss; + struct mc_data cos_factor; struct mc_data shadowed; - double area; - double sun_cos; size_t nb_samples; /* By-receptor data for this entity */ @@ -294,10 +292,8 @@ mc_sampled_init struct mc_sampled* samp) { ASSERT(samp); - samp->cos_loss = MC_DATA_NULL; + samp->cos_factor = MC_DATA_NULL; samp->shadowed = MC_DATA_NULL; - samp->area = 0; - samp->sun_cos = 0; samp->nb_samples = 0; htable_receiver_init(allocator, &samp->mc_rcvs); } @@ -313,10 +309,8 @@ static INLINE res_T mc_sampled_copy(struct mc_sampled* dst, const struct mc_sampled* src) { ASSERT(dst && src); - dst->cos_loss = src->cos_loss; + dst->cos_factor = src->cos_factor; dst->shadowed = src->shadowed; - dst->area = src->area; - dst->sun_cos = src->sun_cos; dst->nb_samples = src->nb_samples; return htable_receiver_copy(&dst->mc_rcvs, &src->mc_rcvs); } @@ -325,10 +319,8 @@ static INLINE res_T mc_sampled_copy_and_release(struct mc_sampled* dst, struct mc_sampled* src) { ASSERT(dst && src); - dst->cos_loss = src->cos_loss; + dst->cos_factor = src->cos_factor; dst->shadowed = src->shadowed; - dst->area = src->area; - dst->sun_cos = src->sun_cos; dst->nb_samples = src->nb_samples; return htable_receiver_copy_and_release(&dst->mc_rcvs, &src->mc_rcvs); } @@ -455,9 +447,12 @@ struct ssol_estimator { size_t failed_count; /* Implicit MC computations */ + struct mc_data cos_factor; + struct mc_data absorbed; struct mc_data shadowed; struct mc_data missing; - struct mc_data cos_loss; /* TODO compute it */ + struct mc_data atmosphere; + struct mc_data reflectivity; struct htable_receiver mc_receivers; /* Per receiver MC */ struct htable_sampled mc_sampled; /* Per sampled instance MC */ @@ -518,21 +513,20 @@ get_mc_sampled struct mc_sampled** out_mc_samp) { struct mc_sampled* mc_samp = NULL; - struct mc_sampled mc_samp_null; res_T res = RES_OK; ASSERT(sampled && inst && out_mc_samp); - mc_sampled_init(inst->dev->allocator, &mc_samp_null); - mc_samp = htable_sampled_find(sampled, &inst); if(!mc_samp) { + struct mc_sampled mc_samp_null; + mc_sampled_init(inst->dev->allocator, &mc_samp_null); res = htable_sampled_set(sampled, &inst, &mc_samp_null); + mc_sampled_release(&mc_samp_null); if(res != RES_OK) goto error; mc_samp = htable_sampled_find(sampled, &inst); } exit: - mc_sampled_release(&mc_samp_null); *out_mc_samp = mc_samp; return res; error: diff --git a/src/ssol_instance.c b/src/ssol_instance.c @@ -203,7 +203,7 @@ ssol_instance_get_area (const struct ssol_instance* instance, double* area) { - if (!instance || !area) return RES_BAD_ARG;; + if (!instance || !area) return RES_BAD_ARG; /* the area of the 3D surface */ *area = instance->shape_rt_area; return RES_OK; diff --git a/src/ssol_mc_receiver.c b/src/ssol_mc_receiver.c @@ -56,9 +56,9 @@ ssol_estimator_get_mc_receiver rcv->Name.SE = rcv->Name.V > 0 ? sqrt(rcv->Name.V / N) : 0; \ } (void)0 SETUP_MC_RESULT(integrated_irradiance); + SETUP_MC_RESULT(integrated_absorbed_irradiance); SETUP_MC_RESULT(absorptivity_loss); SETUP_MC_RESULT(reflectivity_loss); - SETUP_MC_RESULT(cos_loss); #undef SETUP_MC_RESULT rcv->mc__ = mc_rcv1; rcv->N__ = estimator->realisation_count; @@ -106,9 +106,9 @@ ssol_mc_shape_get_mc_primitive prim->Name.SE = 0; \ } (void)0 SETUP_MC_RESULT(integrated_irradiance); + SETUP_MC_RESULT(integrated_absorbed_irradiance); SETUP_MC_RESULT(absorptivity_loss); SETUP_MC_RESULT(reflectivity_loss); - SETUP_MC_RESULT(cos_loss); #undef SETUP_MC_RESULT } else { #define SETUP_MC_RESULT(Name) { \ @@ -119,9 +119,9 @@ ssol_mc_shape_get_mc_primitive prim->Name.SE = prim->Name.V > 0 ? sqrt(prim->Name.V / N) : 0; \ } (void)0 SETUP_MC_RESULT(integrated_irradiance); + SETUP_MC_RESULT(integrated_absorbed_irradiance); SETUP_MC_RESULT(absorptivity_loss); SETUP_MC_RESULT(reflectivity_loss); - SETUP_MC_RESULT(cos_loss); #undef SETUP_MC_RESULT } return RES_OK; diff --git a/src/ssol_object.c b/src/ssol_object.c @@ -23,6 +23,7 @@ #include <rsys/ref_count.h> #include <rsys/rsys.h> #include <rsys/mem_allocator.h> +#include <rsys/double3.h> /******************************************************************************* * Helper functions @@ -139,13 +140,11 @@ ssol_object_add_shaded_shape res = s3d_scene_attach_shape(object->scn_rt, shape->shape_rt); if(res != RES_OK) goto error; mask |= BIT(ATTACH_S3D_RT); - object->scn_rt_area += shape->shape_rt_area; /* Add the shape samp to the sampling scene of the object */ res = s3d_scene_attach_shape(object->scn_samp, shape->shape_samp); if(res != RES_OK) goto error; mask |= BIT(ATTACH_S3D_SAMP); - object->scn_samp_area += shape->shape_samp_area; /* Ask for a shaded shape identifier */ i = darray_shaded_shape_size_get(&object->shaded_shapes); @@ -164,6 +163,9 @@ ssol_object_add_shaded_shape mask |= BIT(REGISTER_SAMP); /* Setup the object shaded shape */ + object->scn_rt_area += shape->shape_rt_area; + object->scn_samp_area += shape->shape_samp_area; + d3_add(object->n, object->n, shape->n); SSOL(shape_ref_get(shape)); SSOL(material_ref_get(front)); SSOL(material_ref_get(back)); @@ -213,6 +215,8 @@ ssol_object_clear(struct ssol_object* obj) htable_shaded_shape_clear(&obj->shaded_shapes_samp); obj->scn_rt_area = 0; + obj->scn_samp_area = 0; + d3_splat(obj->n, 0); S3D(scene_clear(obj->scn_rt)); S3D(scene_clear(obj->scn_samp)); @@ -240,4 +244,3 @@ object_has_shape(struct ssol_object* obj, const struct ssol_shape* shape) S3D(shape_get_id(shape->shape_rt, &id)); return htable_shaded_shape_find(&obj->shaded_shapes_rt, &id) != NULL; } - diff --git a/src/ssol_object_c.h b/src/ssol_object_c.h @@ -48,6 +48,7 @@ struct ssol_object { struct s3d_scene* scn_rt; /* RT scene to instantiate */ struct s3d_scene* scn_samp; /* Sampling scene to instantiate */ double scn_rt_area, scn_samp_area; + double n[3]; /* sum of normals */ struct ssol_device* dev; ref_T ref; diff --git a/src/ssol_shape.c b/src/ssol_shape.c @@ -196,9 +196,10 @@ hyperbol_z (const double p[2], const struct priv_hyperbol_data* hyperbol) { - const double z0 = hyperbol->g_2 + hyperbol->abs_b; + const double z0 = hyperbol->g_square + hyperbol->abs_b; const double r2 = p[0] * p[0] + p[1] * p[1]; - return hyperbol->abs_b * sqrt(1 + r2 * hyperbol->_1_a2) + hyperbol->g_2 - z0; + return hyperbol->abs_b * sqrt(1 + r2 * hyperbol->one_over_a_square) + + hyperbol->g_square - z0; } static FINLINE double @@ -207,7 +208,7 @@ parabol_z const struct priv_parabol_data* parabol) { const double r2 = p[0] * p[0] + p[1] * p[1]; - return r2 * parabol->_1_4f; + return r2 * parabol->one_over_4focal; } static FINLINE double @@ -215,7 +216,7 @@ parabolic_cylinder_z (const double p[2], const struct priv_pcylinder_data* pcyl) { - return (p[1] * p[1]) * pcyl->_1_4f; + return (p[1] * p[1]) * pcyl->one_over_4focal; } static void @@ -474,12 +475,14 @@ mesh_compute_area void (*get_indices)(const unsigned itri, unsigned ids[3], void* data), const unsigned nverts, void (*get_position)(const unsigned ivert, float position[3], void* data), - void* ctx) + void* ctx, + double* normal) { unsigned itri; double area = 0; (void)nverts; + if(normal) d3_splat(normal, 0); FOR_EACH(itri, 0, ntris) { float v0[3], v1[3], v2[3]; double E0[3], E1[3], N[3]; @@ -501,7 +504,9 @@ mesh_compute_area d3_sub(E1, V2, V0); area += d3_len(d3_cross(N, E0, E1)); + if(normal) d3_add(normal, normal, N); } + if (normal) d3_muld(normal, normal, 0.5); return area * 0.5; } @@ -559,7 +564,7 @@ quadric_setup_s3d_shape_rt ASSERT(vdata.get); *rt_area = mesh_compute_area - (ntris, quadric_mesh_get_ids, nverts, vdata.get, &ctx); + (ntris, quadric_mesh_get_ids, nverts, vdata.get, &ctx, NULL); return RES_OK; } @@ -597,7 +602,7 @@ quadric_setup_s3d_shape_samp (shape, ntris, quadric_mesh_get_ids, nverts, &vdata, 1, &ctx); if(res != RES_OK) return res; *samp_area = mesh_compute_area - (ntris, quadric_mesh_get_ids, nverts, quadric_mesh_plane_get_pos, &ctx); + (ntris, quadric_mesh_get_ids, nverts, quadric_mesh_plane_get_pos, &ctx, NULL); return RES_OK; } @@ -709,10 +714,10 @@ quadric_hyperbol_gradient_local { ASSERT(quad && pt && grad); { - const double z0 = quad->g_2 + quad->abs_b; + const double z0 = quad->g_square + quad->abs_b; grad[0] = pt[0]; grad[1] = pt[1]; - grad[2] = -(pt[2] + z0 - quad->g_2) * quad->_a2_b2; + grad[2] = -(pt[2] + z0 - quad->g_square) * quad->a_square_over_b_square; } } @@ -788,14 +793,15 @@ quadric_hyperbol_intersect_local { double dst; const double b2 = quad->abs_b * quad->abs_b; - const double b2_a2 = b2 * quad->_1_a2; - const double z0 = quad->g_2 + quad->abs_b; + const double b2_a2 = b2 * quad->one_over_a_square; + const double z0 = quad->g_square + quad->abs_b; const double a = b2_a2 * (dir[0] * dir[0] + dir[1] * dir[1]) - dir[2] * dir[2]; const double b = - 2 * (b2_a2 * (org[0] * dir[0] + org[1] * dir[1]) - (org[2] + z0 - quad->g_2) * dir[2]); + 2 * (b2_a2 * (org[0] * dir[0] + org[1] * dir[1]) + - (org[2] + z0 - quad->g_square) * dir[2]); const double c = b2_a2 * (org[0] * org[0] + org[1] * org[1]) + b2 - - (org[2] + z0 - quad->g_2) * (org[2] + z0 - quad->g_2); + - (org[2] + z0 - quad->g_square) * (org[2] + z0 - quad->g_square); const int sol = quadric_solve_second(a, b, c, hint, &dst); if(!sol) return 0; @@ -937,7 +943,7 @@ priv_parabol_data_setup { ASSERT(data && parabol); data->focal = parabol->focal; - data->_1_4f = 1 / (4.0 * parabol->focal); + data->one_over_4focal = 1 / (4.0 * parabol->focal); } static FINLINE void @@ -953,10 +959,10 @@ priv_hyperbol_data_setup f = hyperbol->real_focal / g; a2 = g * g * (f - f * f); - data->g_2 = g * 0.5; + data->g_square = g * 0.5; data->abs_b = g * fabs(f - 0.5); - data->_a2_b2 = a2 / (data->abs_b * data->abs_b); - data->_1_a2 = 1 / a2; + data->a_square_over_b_square = a2 / (data->abs_b * data->abs_b); + data->one_over_a_square = 1 / a2; } static FINLINE void @@ -966,7 +972,7 @@ priv_parabolic_cylinder_data_setup { ASSERT(data && parabolic_cylinder); data->focal = parabolic_cylinder->focal; - data->_1_4f = 1 / (4.0 * parabolic_cylinder->focal); + data->one_over_4focal = 1 / (4.0 * parabolic_cylinder->focal); } static INLINE void @@ -1244,6 +1250,7 @@ ssol_punched_surface_setup struct darray_double coords; struct darray_size_t ids; size_t nslices; + double n; res_T res = RES_OK; darray_double_init(shape->dev->allocator, &coords); @@ -1295,6 +1302,10 @@ ssol_punched_surface_setup res = quadric_setup_s3d_shape_samp (psurf->quadric, &coords, &ids, shape->shape_samp, &shape->shape_samp_area); if(res != RES_OK) goto error; + /* the normal to the quadric is known */ + n = psurf->quadric->type == SSOL_QUADRIC_HYPERBOL + ? -shape->shape_samp_area : shape->shape_samp_area; + d3(shape->n, 0, 0, n); exit: darray_double_release(&coords); @@ -1361,7 +1372,7 @@ ssol_mesh_setup (shape->shape_rt, ntris, get_indices, nverts, attrs, nattribs, data); if(res != RES_OK) goto error; shape->shape_rt_area = - mesh_compute_area(ntris, get_indices, nverts, get_position, data); + mesh_compute_area(ntris, get_indices, nverts, get_position, data, shape->n); /* The Star-3D shape to sample is the same of the one to ray-traced */ res = s3d_mesh_copy(shape->shape_rt, shape->shape_samp); diff --git a/src/ssol_shape_c.h b/src/ssol_shape_c.h @@ -28,19 +28,19 @@ enum shape_type { struct priv_parabol_data { double focal; - double _1_4f; + double one_over_4focal; }; struct priv_hyperbol_data { - double g_2; - double _a2_b2; - double _1_a2; + double g_square; + double a_square_over_b_square; + double one_over_a_square; double abs_b; }; struct priv_pcylinder_data { double focal; - double _1_4f; + double one_over_4focal; }; union priv_quadric_data { @@ -57,6 +57,7 @@ struct ssol_shape { union priv_quadric_data priv_quadric; struct ssol_quadric quadric; double shape_rt_area, shape_samp_area; + double n[3]; /* sum of normals */ struct ssol_device* dev; ref_T ref; diff --git a/src/ssol_solver.c b/src/ssol_solver.c @@ -50,9 +50,12 @@ struct thread_context { struct ssp_rng* rng; struct ssf_bsdf* bsdf; + struct mc_data cos_factor; + struct mc_data absorbed; struct mc_data shadowed; struct mc_data missing; - struct mc_data cos_loss; + struct mc_data atmosphere; + struct mc_data reflectivity; struct htable_receiver mc_rcvs; struct htable_sampled mc_samps; struct darray_path paths; /* paths */ @@ -100,9 +103,12 @@ thread_context_copy ASSERT(dst && src); dst->rng = src->rng; dst->bsdf = src->bsdf; + dst->cos_factor = src->cos_factor; + dst->absorbed = src->absorbed; dst->shadowed = src->shadowed; dst->missing = src->missing; - dst->cos_loss = src->cos_loss; + dst->atmosphere = src->atmosphere; + dst->reflectivity = src->reflectivity; res = htable_receiver_copy(&dst->mc_rcvs, &src->mc_rcvs); if(res != RES_OK) return res; res = htable_sampled_copy(&dst->mc_samps, &src->mc_samps); @@ -161,10 +167,12 @@ struct point { double dir[3]; float uv[2]; double wl; /* Sampled wavelength */ - double weight; /* actual weight */ - double absorptivity_loss; - double reflectivity_loss; - double cos_loss; + /* MC weights, before and after hit */ + double incoming_weight, weight; + double cos_factor; /* local cos at the starting point */ + double absorbed_irradiance; /* current hit only */ + double absorptivity_loss_before, absorptivity_loss; + double reflectivity_loss_before, reflectivity_loss; enum ssol_side_flag side; }; @@ -178,7 +186,7 @@ struct point { {0, 0, 0}, /* Direction */ \ {0, 0}, /* UV */ \ 0, /* Wavelength */ \ - 0, 0, 0, 0, /* MC weights */ \ + 0, 0, 0, 0, 0, 0, 0, 0, /* MC weights */ \ SSOL_FRONT /* Side */ \ } static const struct point POINT_NULL = POINT_NULL__; @@ -236,28 +244,25 @@ point_init if(pt->sshape->shape->type != SHAPE_PUNCHED) { 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->cos_factor = surface_sun_cos; } else { - 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, 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); + surface_proxy_cos = fabs(d3_dot(pt->N, tmp_n)); + surface_sun_cos = fabs(d3_dot(tmp_n, pt->dir)); + cos_ratio = 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); + pt->cos_factor = 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->cos_loss.weight += pt->cos_loss; - pt->mc_samp->cos_loss.sqr_weight += pt->cos_loss * pt->cos_loss; pt->mc_samp->nb_samples++; /* Define the primitive side on which the point lies */ @@ -391,13 +396,26 @@ point_shade r = ssf_bsdf_sample(bsdf, rng, wi, frag.Ns, dir, &type, &pdf); ASSERT(0 <= r && r <= 1); } - pt->reflectivity_loss += (1 - r) * pt->weight; - pt->weight *= r; + pt->incoming_weight = pt->weight; + pt->absorptivity_loss_before = pt->absorptivity_loss; + pt->reflectivity_loss_before = pt->reflectivity_loss; + pt->absorbed_irradiance = (1 - r) * pt->weight; + pt->reflectivity_loss += pt->absorbed_irradiance; + pt->weight = pt->incoming_weight - pt->absorbed_irradiance; if(type & SSF_TRANSMISSION) material_get_next_medium(mtl, medium, medium); return RES_OK; } +static FINLINE void +point_hit_virtual(struct point* pt) +{ + pt->absorbed_irradiance = 0; + pt->incoming_weight = pt->weight; + pt->absorptivity_loss_before = pt->absorptivity_loss; + pt->reflectivity_loss_before = pt->reflectivity_loss; +} + static FINLINE int point_is_receiver(const struct point* pt) { @@ -518,9 +536,9 @@ accum_mc_receivers_1side dst->Name.sqr_weight += src->Name.sqr_weight; \ } (void)0 ACCUM_WEIGHT(integrated_irradiance); + ACCUM_WEIGHT(integrated_absorbed_irradiance); ACCUM_WEIGHT(absorptivity_loss); ACCUM_WEIGHT(reflectivity_loss); - ACCUM_WEIGHT(cos_loss); #undef ACCUM_WEIGHT /* Merge the per shape MC */ @@ -555,9 +573,9 @@ accum_mc_receivers_1side mc_prim1_dst->Name.sqr_weight += mc_prim1_src->Name.sqr_weight; \ } (void)0 ACCUM_WEIGHT(integrated_irradiance); + ACCUM_WEIGHT(integrated_absorbed_irradiance); ACCUM_WEIGHT(absorptivity_loss); ACCUM_WEIGHT(reflectivity_loss); - ACCUM_WEIGHT(cos_loss); #undef ACCUM_WEIGHT htable_prim2mc_iterator_next(&it_prim); @@ -585,7 +603,6 @@ accum_mc_sampled(struct mc_sampled* dst, struct mc_sampled* src) dst->Name.weight += src->Name.weight; \ dst->Name.sqr_weight += src->Name.sqr_weight; \ } (void)0 - ACCUM_WEIGHT(cos_loss); ACCUM_WEIGHT(shadowed); #undef ACCUM_WEIGHT @@ -627,43 +644,52 @@ update_mc (const struct point* pt, const size_t irealisation, const size_t ibounce, - struct htable_receiver* mc_rcvs, + struct thread_context* thread_ctx, FILE* output) { struct mc_receiver_1side* mc_rcv1 = NULL; struct mc_receiver_1side* mc_samp_x_rcv1 = NULL; res_T res = RES_OK; - ASSERT(pt && mc_rcvs && point_is_receiver(pt)); - + ASSERT(pt && thread_ctx && point_is_receiver(pt)); + res = point_dump(pt, irealisation, ibounce, output); if(res != RES_OK) goto error; + /* Global MC accumulation */ + #define ACCUM_WEIGHT(Res, W) { \ + Res.weight += (W); \ + Res.sqr_weight += (W)*(W); \ + } (void)0 + ACCUM_WEIGHT(thread_ctx->absorbed, pt->absorbed_irradiance); + #undef ACCUM_WEIGHT + /* Per receiver MC accumulation */ - res = get_mc_receiver_1side(mc_rcvs, pt->inst, pt->side, &mc_rcv1); + res = get_mc_receiver_1side(&thread_ctx->mc_rcvs, pt->inst, pt->side, &mc_rcv1); if(res != RES_OK) goto error; #define ACCUM_WEIGHT(Name, W) { \ mc_rcv1->Name.weight += (W); \ mc_rcv1->Name.sqr_weight += (W)*(W); \ } (void)0 - ACCUM_WEIGHT(integrated_irradiance, pt->weight); - ACCUM_WEIGHT(absorptivity_loss, pt->absorptivity_loss); - ACCUM_WEIGHT(reflectivity_loss, pt->reflectivity_loss); - ACCUM_WEIGHT(cos_loss, pt->cos_loss); + ACCUM_WEIGHT(integrated_irradiance, pt->incoming_weight); + ACCUM_WEIGHT(integrated_absorbed_irradiance, pt->absorbed_irradiance); + ACCUM_WEIGHT(absorptivity_loss, pt->absorptivity_loss_before); + ACCUM_WEIGHT(reflectivity_loss, pt->reflectivity_loss_before); #undef ACCUM_WEIGHT /* Per-sampled/receiver MC accumulation */ res = mc_sampled_get_mc_receiver_1side (pt->mc_samp, pt->inst, pt->side, &mc_samp_x_rcv1); if(res != RES_OK) goto error; + #define ACCUM_WEIGHT(Name, W) { \ mc_samp_x_rcv1->Name.weight += (W); \ mc_samp_x_rcv1->Name.sqr_weight += (W)*(W); \ } (void)0 - ACCUM_WEIGHT(integrated_irradiance, pt->weight); - ACCUM_WEIGHT(absorptivity_loss, pt->absorptivity_loss); - ACCUM_WEIGHT(reflectivity_loss, pt->reflectivity_loss); - ACCUM_WEIGHT(cos_loss, pt->cos_loss); + ACCUM_WEIGHT(integrated_irradiance, pt->incoming_weight); + ACCUM_WEIGHT(integrated_absorbed_irradiance, pt->absorbed_irradiance); + ACCUM_WEIGHT(absorptivity_loss, pt->absorptivity_loss_before); + ACCUM_WEIGHT(reflectivity_loss, pt->reflectivity_loss_before); #undef ACCUM_WEIGHT /* Per primitive receiver MC accumulation */ @@ -677,14 +703,14 @@ 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) { \ - mc_prim1->Name.weight += (W); \ - mc_prim1->Name.sqr_weight += (W)*(W); \ + #define ACCUM_WEIGHT(Name, W) { \ + mc_prim1->Name.weight += (W); \ + mc_prim1->Name.sqr_weight += (W)*(W); \ } (void)0 - ACCUM_WEIGHT(integrated_irradiance, pt->weight); - ACCUM_WEIGHT(absorptivity_loss, pt->absorptivity_loss); - ACCUM_WEIGHT(reflectivity_loss, pt->reflectivity_loss); - ACCUM_WEIGHT(cos_loss, pt->cos_loss); + ACCUM_WEIGHT(integrated_irradiance, pt->incoming_weight); + ACCUM_WEIGHT(integrated_absorbed_irradiance, pt->absorbed_irradiance); + ACCUM_WEIGHT(absorptivity_loss, pt->absorptivity_loss_before); + ACCUM_WEIGHT(reflectivity_loss, pt->reflectivity_loss_before); #undef ACCUM_WEIGHT } @@ -745,16 +771,19 @@ trace_radiative_path 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); + ACCUM_WEIGHT(pt.mc_samp->cos_factor, pt.cos_factor); if(!is_lit) { /* The starting point is not lit */ - pt.mc_samp->shadowed.weight += pt.weight; - pt.mc_samp->shadowed.sqr_weight += pt.weight; - thread_ctx->shadowed.weight += pt.weight; - thread_ctx->shadowed.sqr_weight += pt.weight * pt.weight; + ACCUM_WEIGHT(pt.mc_samp->shadowed, pt.weight); + ACCUM_WEIGHT(thread_ctx->shadowed, pt.weight); if(tracker) path.type = SSOL_PATH_SHADOW; } else { - int hit_a_receiver = 0; - /* Setup the ray as if it starts from the current point position in order * to handle the points that start from a virtual material */ f3_set_d3(org, pt.pos); @@ -765,13 +794,26 @@ trace_radiative_path struct ray_data ray_data = RAY_DATA_NULL; struct ssol_material* mtl; + /* Compute interaction with material */ + mtl = point_get_material(&pt); + if(mtl->type == SSOL_MATERIAL_VIRTUAL) { + point_hit_virtual(&pt); + } else { + /* Modulate the point weight wrt to its scattering functions and + * generate an outgoing direction */ + res = point_shade(&pt, thread_ctx->bsdf, &medium, thread_ctx->rng, pt.dir); + if(res != RES_OK) goto error; + } + if(point_is_receiver(&pt)) { - hit_a_receiver = 1; - res = update_mc(&pt, path_id, depth, &thread_ctx->mc_rcvs, output); + res = update_mc(&pt, path_id, depth, thread_ctx, output); if(res != RES_OK) goto error; } - mtl = point_get_material(&pt); + /* Stop the radiative random walk */ + if(pt.weight == 0) break; + + /* Setup new ray parameters */ if(mtl->type == SSOL_MATERIAL_VIRTUAL) { /* Note that for Virtual materials, the ray parameters 'org' & 'dir' * are not updated to ensure that it pursues its traversal without any @@ -779,15 +821,6 @@ trace_radiative_path range[0] = nextafterf(hit.distance, FLT_MAX); range[1] = FLT_MAX; } else { - /* Modulate the point weight wrt to its scattering functions and - * generate an outgoing direction */ - res = point_shade(&pt, thread_ctx->bsdf, &medium, thread_ctx->rng, pt.dir); - if(res != RES_OK) goto error; - - /* Stop the radiative random walk */ - if(pt.weight == 0) break; - - /* Setup new ray parameters */ f2(range, 0, FLT_MAX); f3_set_d3(org, pt.pos); f3_set_d3(dir, pt.dir); @@ -805,12 +838,12 @@ trace_radiative_path S3D(scene_view_trace_ray(view_rt, org, dir, range, &ray_data, &hit)); if(S3D_HIT_NONE(&hit)) { /* Add the point of the last path segment going to the infinite */ - if(tracker && tracker->infinite_ray_length > 0) { + 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; + if (res != RES_OK) goto error; } break; } @@ -830,15 +863,16 @@ trace_radiative_path if(tracker) { res = path_add_vertex(&path, pt.pos, pt.weight); - if(res != RES_OK) goto error; + 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; - } + ACCUM_WEIGHT(thread_ctx->atmosphere, pt.absorptivity_loss); + /* all the remaining weight is lost */ + ACCUM_WEIGHT(thread_ctx->missing, pt.weight); + #undef ACCUM_WEIGHT + if(tracker) { - path.type = hit_a_receiver ? SSOL_PATH_SUCCESS : SSOL_PATH_MISSING; + path.type = pt.absorbed_irradiance ? SSOL_PATH_SUCCESS : SSOL_PATH_MISSING; } } @@ -846,7 +880,6 @@ trace_radiative_path res = path_register_and_clear(&thread_ctx->paths, &path); if(res != RES_OK) goto error; } - exit: if(tracker) path_release(&path); return res; @@ -967,9 +1000,12 @@ ssol_solve estimator->Name.weight += thread_ctx->Name.weight; \ estimator->Name.sqr_weight += thread_ctx->Name.sqr_weight; \ } (void)0 + ACCUM_WEIGHT(cos_factor); + ACCUM_WEIGHT(absorbed); ACCUM_WEIGHT(shadowed); ACCUM_WEIGHT(missing); - ACCUM_WEIGHT(cos_loss); + ACCUM_WEIGHT(atmosphere); + ACCUM_WEIGHT(reflectivity); #undef ACCUM_WEIGHT } diff --git a/src/test_ssol_instance.c b/src/test_ssol_instance.c @@ -35,7 +35,7 @@ main(int argc, char** argv) struct ssol_vertex_data attrib = SSOL_VERTEX_DATA_NULL; struct ssol_instantiated_shaded_shape sshape; double transform[12] = {1, 0, 0, 0, 1, 0, 0, 0, 1, 10, 0, 0}; - double val[3]; + double val[3], area; size_t n; unsigned i, count; uint32_t id, id1; @@ -79,6 +79,10 @@ main(int argc, char** argv) CHECK(ssol_instance_set_transform(instance, transform), RES_OK); CHECK(ssol_instance_set_transform(instance, transform), RES_OK); + CHECK(ssol_instance_get_area(instance, NULL), RES_BAD_ARG); + CHECK(ssol_instance_get_area(NULL, &area), RES_BAD_ARG); + CHECK(ssol_instance_get_area(instance, &area), RES_OK); + CHECK(ssol_instance_set_receiver(NULL, 0, 0), RES_BAD_ARG); CHECK(ssol_instance_set_receiver(instance, 0, 0), RES_OK); diff --git a/src/test_ssol_object.c b/src/test_ssol_object.c @@ -26,6 +26,7 @@ main(int argc, char** argv) struct ssol_material* mtl; struct ssol_material* mtl2; struct ssol_object* object; + double a; (void) argc, (void) argv; mem_init_proxy_allocator(&allocator, &mem_default_allocator); @@ -73,6 +74,10 @@ main(int argc, char** argv) CHECK(ssol_object_clear(NULL), RES_BAD_ARG); CHECK(ssol_object_clear(object), RES_OK); + CHECK(ssol_object_get_area(object, NULL), RES_BAD_ARG); + CHECK(ssol_object_get_area(NULL, &a), RES_BAD_ARG); + CHECK(ssol_object_get_area(object, &a), RES_OK); + CHECK(ssol_object_ref_put(object), RES_OK); CHECK(ssol_shape_ref_put(shape), RES_OK); CHECK(ssol_shape_ref_put(shape2), RES_OK); diff --git a/src/test_ssol_solver1.c b/src/test_ssol_solver1.c @@ -69,6 +69,7 @@ main(int argc, char** argv) struct ssol_spectrum* abs; struct ssol_atmosphere* atm; struct ssol_estimator* estimator; + struct ssol_mc_sampled sampled; struct ssol_mc_global mc_global; struct ssol_mc_receiver mc_rcv; struct ssol_mc_shape mc_shape; @@ -84,7 +85,7 @@ main(int argc, char** argv) double transform1[12]; /* 3x4 column major matrix */ double transform2[12]; /* 3x4 column major matrix */ double dbl; - size_t i, count, fcount; + size_t i, count, fcount, scount; FILE* tmp = NULL; double m, std; double a_m, a_std; @@ -203,12 +204,12 @@ main(int argc, char** argv) CHECK(ssol_estimator_get_sampled_area(estimator, NULL), RES_BAD_ARG); CHECK(ssol_estimator_get_sampled_area(NULL, &dbl), RES_BAD_ARG); CHECK(ssol_estimator_get_sampled_area(estimator, &dbl), RES_OK); - CHECK(eq_eps(dbl, 12, 1.e-6), 1); + CHECK(eq_eps(dbl, 12, DBL_EPSILON), 1); - CHECK(ssol_estimator_get_count(NULL, NULL), RES_BAD_ARG); - CHECK(ssol_estimator_get_count(estimator, NULL), RES_BAD_ARG); - CHECK(ssol_estimator_get_count(NULL, &count), RES_BAD_ARG); - CHECK(ssol_estimator_get_count(estimator, &count), RES_OK); + CHECK(ssol_estimator_get_realisation_count(NULL, NULL), RES_BAD_ARG); + CHECK(ssol_estimator_get_realisation_count(estimator, NULL), RES_BAD_ARG); + CHECK(ssol_estimator_get_realisation_count(NULL, &count), RES_BAD_ARG); + CHECK(ssol_estimator_get_realisation_count(estimator, &count), RES_OK); CHECK(count, 1); CHECK(ssol_estimator_get_failed_count(NULL, NULL), RES_BAD_ARG); @@ -217,6 +218,17 @@ main(int argc, char** argv) CHECK(ssol_estimator_get_failed_count(estimator, &fcount), RES_OK); CHECK(fcount, 0); + CHECK(ssol_estimator_get_sampled_count(NULL, NULL), RES_BAD_ARG); + CHECK(ssol_estimator_get_sampled_count(estimator, NULL), RES_BAD_ARG); + CHECK(ssol_estimator_get_sampled_count(NULL, &scount), RES_BAD_ARG); + CHECK(ssol_estimator_get_sampled_count(estimator, &scount), RES_OK); + CHECK(scount, 3); + + CHECK(ssol_estimator_get_mc_sampled(NULL, heliostat, &sampled), RES_BAD_ARG); + CHECK(ssol_estimator_get_mc_sampled(estimator, NULL, &sampled), RES_BAD_ARG); + CHECK(ssol_estimator_get_mc_sampled(estimator, heliostat, NULL), RES_BAD_ARG); + CHECK(ssol_estimator_get_mc_sampled(estimator, heliostat, &sampled), RES_OK); + CHECK(ssol_estimator_get_mc_global(NULL, &mc_global), RES_BAD_ARG); CHECK(ssol_estimator_get_mc_global(estimator, NULL), RES_BAD_ARG); CHECK(ssol_estimator_get_mc_global(estimator, &mc_global), RES_OK); @@ -232,6 +244,7 @@ main(int argc, char** argv) CHECK(ssol_instance_sample(secondary, 0), RES_OK); CHECK(ssol_instance_sample(heliostat, 0), RES_OK); 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); CHECK(ssol_instance_sample(secondary, 1), RES_OK); @@ -290,7 +303,7 @@ main(int argc, char** argv) #define GET_MC_GLOBAL ssol_estimator_get_mc_global 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(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); CHECK(fclose(tmp), 0); @@ -308,8 +321,9 @@ main(int argc, char** argv) CHECK(ssol_estimator_get_mc_global(estimator, &mc_global), RES_OK); printf("Shadows = %g +/- %g; ", mc_global.shadowed.E, mc_global.shadowed.SE); printf("Missing = %g +/- %g; ", mc_global.missing.E, mc_global.missing.SE); + printf("Cos = %g +/- %g; ", mc_global.cos_factor.E, mc_global.cos_factor.SE); CHECK(eq_eps(mc_global.shadowed.E, m, 2 * dbl), 1); - CHECK(eq_eps(mc_global.missing.E, m, 2*mc_global.missing.SE), 1); + CHECK(eq_eps(mc_global.missing.E, 2*m, 2*mc_global.missing.SE), 1); CHECK(GET_MC_RCV(NULL, NULL, SSOL_BACK, NULL), RES_BAD_ARG); CHECK(GET_MC_RCV(estimator, NULL, SSOL_BACK, NULL), RES_BAD_ARG); CHECK(GET_MC_RCV(NULL, target, SSOL_BACK, NULL), RES_BAD_ARG); @@ -338,7 +352,7 @@ main(int argc, char** argv) NCHECK(tmp = tmpfile(), 0); CHECK(ssol_solve(scene, rng, N__, 0, tmp, &estimator), RES_OK); - CHECK(ssol_estimator_get_count(estimator, &count), 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); CHECK(fclose(tmp), 0); @@ -348,8 +362,10 @@ main(int argc, char** argv) CHECK(ssol_estimator_get_mc_global(estimator, &mc_global), RES_OK); printf("Shadows = %g +/- %g; ", mc_global.shadowed.E, mc_global.shadowed.SE); printf("Missing = %g +/- %g; ", mc_global.missing.E, mc_global.missing.SE); + printf("Cos = %g +/- %g; ", mc_global.cos_factor.E, mc_global.cos_factor.SE); CHECK(eq_eps(mc_global.shadowed.E, 0, 1e-4), 1); - CHECK(eq_eps(mc_global.missing.E, 0, 1e-4), 1); + CHECK(eq_eps(mc_global.missing.E, m, 1e-4), 1); + CHECK(eq_eps(mc_global.cos_factor.E, COS, 1e-4), 1); CHECK(GET_MC_RCV(estimator, target, SSOL_FRONT, &mc_rcv), RES_OK); printf("Ir(target) = %g +/- %g\n", mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE); @@ -369,7 +385,7 @@ main(int argc, char** argv) NCHECK(tmp = tmpfile(), 0); CHECK(ssol_solve(scene, rng, N__, 0, tmp, &estimator), RES_OK); - CHECK(ssol_estimator_get_count(estimator, &count), 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); CHECK(fclose(tmp), 0); @@ -382,13 +398,16 @@ main(int argc, char** argv) CHECK(ssol_estimator_get_mc_global(estimator, &mc_global), RES_OK); printf("Shadows = %g +/- %g; ", mc_global.shadowed.E, mc_global.shadowed.SE); printf("Missing = %g +/- %g; ", mc_global.missing.E, mc_global.missing.SE); + printf("Cos = %g +/- %g; ", mc_global.cos_factor.E, mc_global.cos_factor.SE); CHECK(eq_eps(mc_global.shadowed.E, 0, 1e-4), 1); - CHECK(eq_eps(mc_global.missing.E, 0, 1e-4), 1); + CHECK(eq_eps(mc_global.missing.E, m, 1e-4), 1); + CHECK(eq_eps(mc_global.cos_factor.E, COS, 1e-4), 1); CHECK(GET_MC_RCV(estimator, 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, m, 1e-8), 1); CHECK(eq_eps(mc_rcv.integrated_irradiance.SE, std, 1e-4), 1); + CHECK(eq_eps(mc_global.cos_factor.E, COS, 1e-4), 1); CHECK(ssol_estimator_ref_put(estimator), RES_OK); /* Check atmosphere model and imperfect mirror: there are losses */ @@ -417,7 +436,7 @@ main(int argc, char** argv) NCHECK(tmp = tmpfile(), 0); CHECK(ssol_solve(scene, rng, N__, 0, tmp, &estimator), RES_OK); - CHECK(ssol_estimator_get_count(estimator, &count), 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); CHECK(fclose(tmp), 0); @@ -428,9 +447,13 @@ main(int argc, char** argv) CHECK(eq_eps(a_std, 0, 1e-4), 1); CHECK(ssol_estimator_get_mc_global(estimator, &mc_global), RES_OK); printf("Shadows = %g +/- %g; ", mc_global.shadowed.E, mc_global.shadowed.SE); - printf("Missing = %g +/- %g\n", mc_global.missing.E, mc_global.missing.SE); + printf("Missing = %g +/- %g; ", mc_global.missing.E, mc_global.missing.SE); + printf("Atmosphere = %g +/- %g; ", mc_global.atmosphere.E, mc_global.atmosphere.SE); + printf("Cos = %g +/- %g\n", mc_global.cos_factor.E, mc_global.cos_factor.SE); CHECK(eq_eps(mc_global.shadowed.E, 0, 1e-4), 1); - CHECK(eq_eps(mc_global.missing.E, 0, 1e-4), 1); + CHECK(eq_eps(mc_global.missing.E + mc_global.atmosphere.E + mc_global.absorbed.E, + m, 1e-4), 1); + CHECK(eq_eps(mc_global.cos_factor.E, COS, 1e-4), 1); CHECK(GET_MC_RCV(estimator, target, SSOL_FRONT, &mc_rcv), RES_OK); printf ("\tIr(target) = %g +/- %g (%.2g %%)\n", @@ -447,11 +470,11 @@ main(int argc, char** argv) mc_rcv.reflectivity_loss.E, mc_rcv.reflectivity_loss.SE, 100 * mc_rcv.reflectivity_loss.E / m); - printf - ("\tCos Loss(target) = %g +/- %g (%.2g %%)\n", - mc_rcv.cos_loss.E, - mc_rcv.cos_loss.SE, - 100 * mc_rcv.cos_loss.E / m); + + CHECK(ssol_estimator_get_sampled_count(estimator, &scount), RES_OK); + CHECK(ssol_estimator_get_mc_sampled(estimator, heliostat, &sampled), RES_BAD_ARG); + CHECK(ssol_estimator_get_mc_sampled(estimator, heliostat2, &sampled), RES_OK); + CHECK(eq_eps(mc_rcv.integrated_irradiance.E, a_m, 1e-8), 1); CHECK(eq_eps(mc_rcv.integrated_irradiance.SE, a_std, 1e-4), 1); CHECK(eq_eps @@ -462,9 +485,7 @@ main(int argc, char** argv) ( mc_rcv.integrated_irradiance.E + mc_rcv.absorptivity_loss.E + mc_rcv.reflectivity_loss.E - + mc_rcv.cos_loss.E, 4 * DNI, 1e-8), 1); - CHECK(eq_eps(mc_rcv.cos_loss.E / (4 * DNI), 1 - COS, 1e-8), 1); - + + (1 - mc_global.cos_factor.E) * 4 * DNI, 4 * DNI, 1e-8), 1); CHECK(ssol_mc_receiver_get_mc_shape(NULL, NULL, NULL), RES_BAD_ARG); CHECK(ssol_mc_receiver_get_mc_shape(&mc_rcv, NULL, NULL), RES_BAD_ARG); CHECK(ssol_mc_receiver_get_mc_shape(NULL, square, NULL), RES_BAD_ARG); @@ -517,7 +538,7 @@ main(int argc, char** argv) CHECK(ssol_spectrum_setup(abs, get_wlen, 2, &desc), RES_OK); NCHECK(tmp = tmpfile(), 0); CHECK(ssol_solve(scene, rng, N__, 0, tmp, &estimator), RES_OK); - CHECK(ssol_estimator_get_count(estimator, &count), 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); CHECK(fclose(tmp), 0); @@ -528,8 +549,10 @@ main(int argc, char** argv) CHECK(ssol_estimator_get_mc_global(estimator, &mc_global), RES_OK); printf("Shadows = %g +/- %g; ", mc_global.shadowed.E, mc_global.shadowed.SE); printf("Missing = %g +/- %g; ", mc_global.missing.E, mc_global.missing.SE); + printf("Cos = %g +/- %g; ", mc_global.cos_factor.E, mc_global.cos_factor.SE); CHECK(eq_eps(mc_global.shadowed.E, 0, 1e-4), 1); - CHECK(eq_eps(mc_global.missing.E, 0, 1e-4), 1); + CHECK(eq_eps(mc_global.missing.E, m, 1e-4), 1); + CHECK(eq_eps(mc_global.cos_factor.E, COS, 1e-4), 1); CHECK(GET_MC_RCV(estimator, 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_solver2.c b/src/test_ssol_solver2.c @@ -182,20 +182,23 @@ main(int argc, char** argv) #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(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); CHECK(fclose(tmp), 0); printf("Ir = %g +/- %g\n", m, std); -#define DNI_cos (1000 * cos(PI / 4)) +#define COS cos(PI / 4) +#define DNI_cos (1000 * COS) CHECK(eq_eps(m, 4 * DNI_cos, 4 * DNI_cos * 1e-4), 1); #define SQR(x) ((x)*(x)) CHECK(eq_eps(std, 0, 1e-4), 1); CHECK(ssol_estimator_get_mc_global(estimator, &mc_global), RES_OK); printf("Shadows = %g +/- %g\n", mc_global.shadowed.E, mc_global.shadowed.SE); - printf("Missing = %g +/- %g", mc_global.missing.E, mc_global.missing.SE); + printf("Missing = %g +/- %g\n", mc_global.missing.E, mc_global.missing.SE); + printf("Cos = %g +/- %g\n", mc_global.cos_factor.E, mc_global.cos_factor.SE); CHECK(eq_eps(mc_global.shadowed.E, 0, 1e-4), 1); - CHECK(eq_eps(mc_global.missing.E, 0, 1e-4), 1); + CHECK(eq_eps(mc_global.missing.E, 4 * DNI_cos, 1e-4), 1); /* nothing absorbed */ + CHECK(eq_eps(mc_global.cos_factor.E, COS, 1e-4), 1); CHECK(GET_MC_RCV(estimator, heliostat1, SSOL_BACK, &mc_rcv), RES_BAD_ARG); CHECK(GET_MC_RCV(estimator, secondary, SSOL_FRONT, &mc_rcv), RES_OK); printf("Ir(secondary) = %g +/- %g\n", diff --git a/src/test_ssol_solver2b.c b/src/test_ssol_solver2b.c @@ -185,12 +185,13 @@ main(int argc, char** argv) #define N__ 50000 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(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); CHECK(fclose(tmp), 0); printf("Ir = %g +/- %g\n", m, std); -#define DNI_cos (1000 * cos(PI / 4)) +#define COS cos(PI / 4) +#define DNI_cos (1000 * COS) CHECK(eq_eps(m, 2 * DNI_cos, MMAX(2 * DNI_cos * 1e-2, std)), 1); #define SQR(x) ((x)*(x)) CHECK(eq_eps(std, @@ -198,8 +199,10 @@ main(int argc, char** argv) CHECK(ssol_estimator_get_mc_global(estimator, &mc_global), RES_OK); 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); + printf("Cos = %g +/- %g\n", mc_global.cos_factor.E, mc_global.cos_factor.SE); CHECK(eq_eps(mc_global.shadowed.E, 0, 1e-4), 1); - CHECK(eq_eps(mc_global.missing.E, 0, 1e-4), 1); + CHECK(eq_eps(mc_global.missing.E, 4 * DNI_cos, 1e-4), 1); /* nothing absorbed */ + CHECK(eq_eps(mc_global.cos_factor.E, COS, 1e-4), 1); CHECK(ssol_estimator_get_mc_receiver (estimator, target, SSOL_FRONT, &mc_rcv), RES_OK); printf("Ir(target) = %g +/- %g\n", diff --git a/src/test_ssol_solver3.c b/src/test_ssol_solver3.c @@ -70,6 +70,7 @@ main(int argc, char** argv) struct ssol_mc_receiver mc_rcv; double dir[3]; double transform[12]; /* 3x4 column major matrix */ + double area; size_t count; FILE* tmp; double m, std; @@ -139,12 +140,13 @@ main(int argc, char** argv) #define N__ 20000 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(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); CHECK(fclose(tmp), 0); printf("Ir = %g +/- %g\n", m, std); -#define DNI_cos (1000 * cos(PI / 4)) +#define COS cos(PI / 4) +#define DNI_cos (1000 * COS) CHECK(eq_eps(m, 4 * DNI_cos, 4 * DNI_cos * 2e-1), 1); #define SQR(x) ((x)*(x)) CHECK(eq_eps(std, @@ -152,8 +154,10 @@ main(int argc, char** argv) CHECK(ssol_estimator_get_mc_global(estimator, &mc_global), RES_OK); 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); + printf("Cos = %g +/- %g\n", mc_global.cos_factor.E, mc_global.cos_factor.SE); CHECK(eq_eps(mc_global.shadowed.E, 0, 1e-4), 1); - CHECK(eq_eps(mc_global.missing.E, 0, 1e-4), 1); + CHECK(eq_eps(mc_global.missing.E, 400 * DNI_cos, 1e-4), 1); /* nothing absorbed */ + CHECK(eq_eps(mc_global.cos_factor.E, COS, 1e-4), 1); CHECK(ssol_estimator_get_mc_receiver (estimator, target, SSOL_FRONT, &mc_rcv), RES_OK); printf("Ir(target) = %g +/- %g\n", @@ -162,6 +166,8 @@ main(int argc, char** argv) CHECK(eq_eps(mc_rcv.integrated_irradiance.SE, std, 1e-4), 1); CHECK(ssol_estimator_get_failed_count(estimator, &count), RES_OK); CHECK(count, 0); + CHECK(ssol_instance_get_area(heliostat, &area), RES_OK); + CHECK(eq_eps(area, 400, DBL_EPSILON), 1); /* Free data */ CHECK(ssol_instance_ref_put(heliostat), RES_OK); diff --git a/src/test_ssol_solver4.c b/src/test_ssol_solver4.c @@ -112,6 +112,7 @@ main(int argc, char** argv) carving.context = &POLY_EDGES__; quadric.type = SSOL_QUADRIC_PARABOL; quadric.data.parabol.focal = FOCAL; + quadric.slices_count_hint = 100; punched.nb_carvings = 1; punched.quadric = &quadric; punched.carvings = &carving; @@ -144,18 +145,19 @@ main(int argc, char** argv) CHECK(ssol_scene_attach_instance(scene, target2), RES_OK); NCHECK(tmp = tmpfile(), 0); -#define N__ 10000 +#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(target1, &r_id1), RES_OK); CHECK(ssol_instance_get_id(target2, &r_id2), RES_OK); - CHECK(ssol_estimator_get_count(estimator, &count), RES_OK); + CHECK(ssol_estimator_get_realisation_count(estimator, &count), RES_OK); CHECK(count, N__); CHECK(pp_sum(tmp, (int32_t)r_id1, count, &m1, &std1), RES_OK); CHECK(pp_sum(tmp, (int32_t)r_id2, count, &m2, &std2), RES_OK); CHECK(fclose(tmp), 0); printf("Ir = %g +/- %g\n", m1, std1); -#define DNI_cos (1000 * cos(0)) +#define COS cos(0) +#define DNI_cos (1000 * COS) CHECK(eq_eps(m1, 400 * DNI_cos, 400 * DNI_cos * 1e-4), 1); CHECK(eq_eps(std1, 0, 1), 1); CHECK(m1, m2); @@ -163,18 +165,19 @@ main(int argc, char** argv) CHECK(ssol_estimator_get_mc_global(estimator, &mc_global), RES_OK); 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); + printf("Cos = %g +/- %g\n", mc_global.cos_factor.E, mc_global.cos_factor.SE); CHECK(eq_eps(mc_global.shadowed.E, 0, 1e-4), 1); - CHECK(eq_eps(mc_global.missing.E, 0, 1e-4), 1); + CHECK(eq_eps(mc_global.missing.E, 400 * DNI_cos, 1e-2), 1); /* nothing absorbed */ CHECK(GET_MC_RCV(estimator, target1, 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, m1, 1e-8), 1); - CHECK(eq_eps(mc_rcv.integrated_irradiance.SE, std1, 1e-4), 1); + CHECK(eq_eps(mc_rcv.integrated_irradiance.E, m1, 1e-2), 1); + CHECK(eq_eps(mc_rcv.integrated_irradiance.SE, std1, 1e-2), 1); CHECK(GET_MC_RCV(estimator, target2, SSOL_FRONT, &mc_rcv), RES_OK); printf("Ir(target2) = %g +/- %g\n", mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE); - CHECK(eq_eps(mc_rcv.integrated_irradiance.E, m2, 1e-8), 1); - CHECK(eq_eps(mc_rcv.integrated_irradiance.SE, std2, 1e-4), 1); + CHECK(eq_eps(mc_rcv.integrated_irradiance.E, m2, 1e-2), 1); + CHECK(eq_eps(mc_rcv.integrated_irradiance.SE, std2, 1e-2), 1); CHECK(ssol_estimator_get_failed_count(estimator, &count), RES_OK); CHECK(count, 0); diff --git a/src/test_ssol_solver5.c b/src/test_ssol_solver5.c @@ -140,19 +140,21 @@ main(int argc, char** argv) #define N__ 10000 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(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); CHECK(fclose(tmp), 0); printf("Ir = %g +/- %g\n", m, std); -#define DNI_cos (1000 * cos(0)) +#define COS cos(0) +#define DNI_cos (1000 * COS) CHECK(eq_eps(m, 400 * DNI_cos, 20), 1); CHECK(eq_eps(std, 0, 1), 1); CHECK(ssol_estimator_get_mc_global(estimator, &mc_global), RES_OK); 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); + printf("Cos = %g +/- %g\n", mc_global.cos_factor.E, mc_global.cos_factor.SE); CHECK(eq_eps(mc_global.shadowed.E, 0, 1e-4), 1); - CHECK(eq_eps(mc_global.missing.E, 0, 1e-4), 1); + CHECK(eq_eps(mc_global.missing.E, 400 * DNI_cos, 1e-4), 1); /* nothing absorbed */ CHECK(ssol_estimator_get_mc_receiver (estimator, target, SSOL_FRONT, &mc_rcv), RES_OK); printf("Ir(target) = %g +/- %g\n", diff --git a/src/test_ssol_solver6.c b/src/test_ssol_solver6.c @@ -181,8 +181,10 @@ main(int argc, char** argv) 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("Absorbed = %g +/- %g\n", mc_global.absorbed.E, mc_global.absorbed.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); + printf("Cos = %g +/- %g\n", mc_global.cos_factor.E, mc_global.cos_factor.SE); CHECK(eq_eps(mc_global.shadowed.E, 100000, 2 * 100000/sqrt(N__)), 1); CHECK(eq_eps(mc_global.missing.E, 0, 0), 1); @@ -190,11 +192,13 @@ main(int argc, char** argv) printf("Ir(target1) = %g +/- %g\n", mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE); CHECK(eq_eps(mc_rcv.integrated_irradiance.E, 100000, 2*100000/sqrt(N__)), 1); + CHECK(mc_rcv.integrated_irradiance.E, mc_rcv.integrated_absorbed_irradiance.E); CHECK(GET_MC_RCV(estimator, target2, SSOL_BACK, &mc_rcv), RES_OK); printf("Ir(target2) = %g +/- %g\n", mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE); CHECK(eq_eps(mc_rcv.integrated_irradiance.E, 0, 1), 1); + CHECK(mc_rcv.integrated_irradiance.E, mc_rcv.integrated_absorbed_irradiance.E); /* Free data */ CHECK(ssol_instance_ref_put(heliostat1), RES_OK); diff --git a/src/test_ssol_solver7.c b/src/test_ssol_solver7.c @@ -19,12 +19,7 @@ #include <rsys/mem_allocator.h> #include <rsys/image.h> -#define SCREEN_GAMMA 2.2 -#define WIDTH 640 -#define HEIGHT 480 -#define PROJ_RATIO (double)WIDTH/(double)HEIGHT - -#define REFLECTIVITY 0 +#define REFLECTIVITY 0.1 #include "test_ssol_materials.h" #define TARGET_SZ 2 @@ -85,7 +80,7 @@ main(int argc, char** argv) struct ssol_estimator* estimator; struct ssol_mc_global mc_global; struct ssol_mc_receiver mc_rcv; - double dir[3]; + double dir[3], area; double transform[12]; /* 3x4 column major matrix */ FILE* tmp; /* primary is a parabol */ @@ -201,20 +196,26 @@ main(int argc, char** argv) CHECK(ssol_solve(scene, rng, N__, 0, tmp, &estimator), RES_OK); CHECK(fclose(tmp), 0); - printf("Total = %g\n", TOTAL); CHECK(ssol_estimator_get_mc_global(estimator, &mc_global), RES_OK); - printf("Shadows = %g +/- %g\n", - mc_global.shadowed.E, mc_global.shadowed.SE); + CHECK(ssol_estimator_get_sampled_area(estimator, &area), RES_OK); + printf("Total = %g\n", area * DNI_cos); + CHECK(eq_eps(area * DNI_cos, TOTAL, TOTAL * 1e-4), 1); + printf("Absorbed = %g +/- %g\n", mc_global.absorbed.E, mc_global.absorbed.SE); + printf("Cos = %g +/- %g\n", mc_global.cos_factor.E, mc_global.cos_factor.SE); + printf("Shadows = %g +/- %g\n", mc_global.shadowed.E, mc_global.shadowed.SE); CHECK(eq_eps(mc_global.shadowed.E, 0, 1e-4), 1); - printf("Missing = %g +/- %g\n", - mc_global.missing.E, mc_global.missing.SE); - CHECK(eq_eps(mc_global.missing.E, 0, 1e-4), 1); + printf("Missing = %g +/- %g\n", mc_global.missing.E, mc_global.missing.SE); CHECK(GET_MC_RCV(estimator, target, SSOL_FRONT, &mc_rcv), RES_OK); + printf("Abs(target1) = %g +/- %g\n", + mc_rcv.integrated_absorbed_irradiance.E, + mc_rcv.integrated_absorbed_irradiance.SE); printf("Ir(target1) = %g +/- %g\n", mc_rcv.integrated_irradiance.E, mc_rcv.integrated_irradiance.SE); CHECK(eq_eps(mc_rcv.integrated_irradiance.E, TOTAL, TOTAL * 1e-4), 1); - CHECK(eq_eps(mc_rcv.integrated_irradiance.SE, 0, 1e-4), 1); + CHECK(eq_eps(mc_rcv.integrated_absorbed_irradiance.E, + (1 - REFLECTIVITY) * TOTAL, (1 - REFLECTIVITY) *TOTAL * 1e-4), 1); + CHECK(eq_eps(mc_rcv.integrated_irradiance.SE, 0, 1e-2), 1); /* Free data */ CHECK(ssol_instance_ref_put(heliostat), RES_OK); diff --git a/src/test_ssol_solver8.c b/src/test_ssol_solver8.c @@ -144,19 +144,19 @@ main(int argc, char** argv) #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(ssol_estimator_get_realisation_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("Cos = %g +/- %g\n", mc_global.cos_factor.E, mc_global.cos_factor.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(eq_eps(mc_global.missing.E, S * DNI, + 3 * mc_global.missing.SE), 1); /* nothing absorbed */ 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); diff --git a/src/test_ssol_solver9.c b/src/test_ssol_solver9.c @@ -138,7 +138,7 @@ main(int argc, char** argv) #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(ssol_estimator_get_realisation_count(estimator, &count), RES_OK); CHECK(count, N__); CHECK(fclose(tmp), 0); CHECK(ssol_estimator_get_failed_count(estimator, &count), RES_OK); @@ -146,13 +146,11 @@ main(int argc, char** argv) #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("Cos = %g +/- %g\n", mc_global.cos_factor.E, mc_global.cos_factor.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), + CHECK(eq_eps(mc_global.missing.E, MMAX(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",