solstice-solver

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

commit 718eab90f6d20ac85bc702ab921f708ed60bd311
parent ee44c5cd546e95227cc9d66afd9f19379295d355
Author: Christophe Coustet <christophe.coustet@meso-star.com>
Date:   Fri,  7 Oct 2016 18:39:13 +0200

Add global indicators (always-ON MC computations)

The success_mask concept is removed: success = some receiver has been hit

Diffstat:
Mcmake/CMakeLists.txt | 3+++
Msrc/ssol.h | 45++++++++++++++++++++++++++++++++++++++++++++-
Msrc/ssol_scene.c | 7-------
Msrc/ssol_solver.c | 88++++++++++++++++++++++++++++++++++++++++----------------------------------------
Msrc/ssol_solver_c.h | 25++-----------------------
Msrc/test_ssol_solver1.c | 93+++++++++++++++++++++++++++++++++++++++++++++++++++++++++----------------------
Msrc/test_ssol_solver2.c | 18++++++++++++++----
Msrc/test_ssol_solver2b.c | 20+++++++++++++++-----
Msrc/test_ssol_solver3.c | 17++++++++++++++---
Msrc/test_ssol_solver3N.c | 15+++++++++++++--
Msrc/test_ssol_solver4.c | 17++++++++++++++---
Msrc/test_ssol_solver5.c | 15+++++++++++++--
12 files changed, 243 insertions(+), 120 deletions(-)

diff --git a/cmake/CMakeLists.txt b/cmake/CMakeLists.txt @@ -55,6 +55,7 @@ set(VERSION ${VERSION_MAJOR}.${VERSION_MINOR}.${VERSION_PATCH}) set(SSOL_FILES_SRC ssol_atmosphere.c ssol_device.c + ssol_estimator.c ssol_image.c ssol_material.c ssol_object.c @@ -74,6 +75,7 @@ set(SSOL_FILES_INC ssol_atmosphere_c.h ssol_c.h ssol_device_c.h + ssol_estimator_c.h ssol_image_c.h ssol_material_c.h ssol_object_c.h @@ -131,6 +133,7 @@ if(NOT NO_TEST) new_test(test_ssol_atmosphere) new_test(test_ssol_device) + new_test(test_ssol_estimator) new_test(test_ssol_image) new_test(test_ssol_material) new_test(test_ssol_object) diff --git a/src/ssol.h b/src/ssol.h @@ -56,6 +56,7 @@ struct ssol_scene; struct ssol_shape; struct ssol_spectrum; struct ssol_sun; +struct ssol_estimator; enum ssol_clipping_op { SSOL_AND, @@ -209,6 +210,21 @@ struct receiver_data { float uv[2]; }; +/* result for MC simulations */ +struct ssol_estimator_status { + double E; /* Expected value */ + double V; /* Variance */ + double SE; /* Standard error, i.e. sqrt(V / N) */ + size_t N; /* Samples count */ +}; + +/* the always-ON indicators (MC computations) */ +enum status_type { + STATUS_SHADOW, + STATUS_MISSING, + STATUS_TYPES_COUNT__ +}; + /* * All the ssol structures are ref counted. Once created with the appropriated * `ssol_<TYPE>_create' function, the caller implicitly owns the created data, @@ -571,6 +587,32 @@ ssol_atmosphere_set_uniform_absorbtion struct ssol_spectrum* spectrum); /******************************************************************************* +* Estimator API - Describe the state of a simulation. +******************************************************************************/ +SSOL_API res_T +ssol_estimator_create + (struct ssol_device* dev, + struct ssol_estimator** estimator); + +SSOL_API res_T +ssol_estimator_get_status + (const struct ssol_estimator* estimator, + enum status_type type, + struct ssol_estimator_status* status); + +SSOL_API res_T +ssol_estimator_clear + (struct ssol_estimator* estimator); + +SSOL_API res_T +ssol_estimator_ref_get + (struct ssol_estimator* estimator); + +SSOL_API res_T +ssol_estimator_ref_put + (struct ssol_estimator* estimator); + +/******************************************************************************* * Miscellaneous functions ******************************************************************************/ SSOL_API res_T @@ -578,7 +620,8 @@ ssol_solve (struct ssol_scene* scn, struct ssp_rng* rng, const size_t realisations_count, - FILE* output); + FILE* output, + struct ssol_estimator* estimator); END_DECLS diff --git a/src/ssol_scene.c b/src/ssol_scene.c @@ -434,13 +434,6 @@ hit_filter_function &rs->data.receiver_record_candidates, &candidate); } - /* register success mask */ - if(seg->hit_front) { - rs->success_mask |= inst->target_front_mask; - } - else { - rs->success_mask |= inst->target_back_mask; - } if(seg->hit_material->type == MATERIAL_VIRTUAL) { return 1; /* Discard virtual material */ } diff --git a/src/ssol_solver.c b/src/ssol_solver.c @@ -18,6 +18,7 @@ #include "ssol_atmosphere_c.h" #include "ssol_solver_c.h" #include "ssol_device_c.h" +#include "ssol_estimator_c.h" #include "ssol_scene_c.h" #include "ssol_shape_c.h" #include "ssol_object_c.h" @@ -38,11 +39,6 @@ #include <star/ssf.h> -#define END_TEXT__ \ - { "NONE", "SUCCESS", "SHADOW", "POINTING", "MISSING", "BLOCKED", "ERROR" } - -static const char* END_TEXT[] = END_TEXT__; - /******************************************************************************* * Helper functions ******************************************************************************/ @@ -289,7 +285,7 @@ setup_next_segment(struct realisation* rs) if (++rs->s_idx >= darray_segment_size_get(&rs->segments)) { res = darray_segment_resize(&rs->segments, rs->s_idx + 1); - if (res != RES_OK) return res; + if (res != RES_OK) goto error; } prev = previous_segment(rs); seg = current_segment(rs); @@ -309,10 +305,7 @@ setup_next_segment(struct realisation* rs) res = material_shade (prev->hit_material, &data->fragment, rs->wavelength, data->bsdf); - if (res != RES_OK) { - rs->end = TERM_ERR; - return res; - } + if (res != RES_OK) goto error; /* By convention, Star-SF assumes that incoming and reflected directions * point outward the surface => negate incoming dir */ @@ -328,7 +321,13 @@ setup_next_segment(struct realisation* rs) (rs->data.scene->atmosphere, prev->hit_distance, rs->wavelength); } + end: return res; + + error: + rs->end = 1; + rs->error = 1; + goto end; } void @@ -398,10 +397,11 @@ static void reset_realisation(size_t cpt, struct realisation* rs) { rs->s_idx = 0; - rs->end = TERM_NONE; - rs->mode = MODE_STD; + rs->end = 0; + rs->error = 0; + rs->shadow = 0; + rs->success = 0; rs->rs_id = cpt; - rs->success_mask = 0; reset_starting_point(&rs->start); SSF(bsdf_clear(rs->data.bsdf)); /* reset first segment (always used) */ @@ -596,12 +596,14 @@ receive_sunlight(struct realisation* rs) seg->self_instance = start->instance; seg->self_front = start->front_exposed; seg->sun_segment = 1; + seg->weight = rs->data.sampled_area * sun->dni * start->cos_sun; /* search for occlusions from starting point */ ASSERT(rs->s_idx == 0); /* sun segment */ solstice_trace_ray(rs); if (!S3D_HIT_NONE(&seg->hit)) { - rs->end = TERM_SHADOW; + rs->end = 1; + rs->shadow = 1; return 0; } @@ -620,7 +622,6 @@ receive_sunlight(struct realisation* rs) seg->self_instance = NULL; seg->hit_front = seg->self_front; seg->self_front = NON_BOOL; - seg->weight = rs->data.sampled_area * sun->dni * start->cos_sun; ASSERT(seg->weight > 0); /* fill fragment from starting point */ @@ -654,14 +655,6 @@ receive_sunlight(struct realisation* rs) &rs->data.receiver_record_candidates, &candidate); } - /* register success mask (normal orientation has already been checked) */ - if (start->front_exposed) { - rs->success_mask |= start->instance->target_front_mask; - } - else { - rs->success_mask |= start->instance->target_back_mask; - } - return 1; } @@ -732,6 +725,7 @@ filter_receiver_hit_candidates(struct realisation* rs) out.weight = weight; f2_set(out.uv, candidates->uv); fwrite(&out, sizeof(struct receiver_data), 1, rs->data.out_stream); + rs->success = 1; } /* reset candidates */ darray_receiver_record_clear(&rs->data.receiver_record_candidates); @@ -752,7 +746,7 @@ propagate(struct realisation* rs) filter_receiver_hit_candidates(rs); if (S3D_HIT_NONE(&seg->hit)) { - rs->end = TERM_MISSING; + rs->end = 1; return; } check_segment(seg); @@ -770,13 +764,14 @@ ssol_solve (struct ssol_scene* scene, struct ssp_rng* rng, const size_t realisations_count, - FILE* output) + FILE* output, + struct ssol_estimator* estimator) { struct realisation rs; size_t r; res_T res = RES_OK; - if (!scene || !rng || !output || !realisations_count) + if (!scene || !rng || !output || !estimator || !realisations_count) return RES_BAD_ARG; res = check_scene(scene); @@ -788,6 +783,8 @@ ssol_solve for (r = 0; r < realisations_count; ) { do { + struct segment* seg; + double w; reset_realisation(r, &rs); sample_starting_point(&rs); sample_input_sundir(&rs); @@ -797,32 +794,35 @@ ssol_solve if (receive_sunlight(&rs)) { /* start propagating from mirror */ do { - if (RES_OK != setup_next_segment(&rs)) { - rs.end = TERM_ERR; - } - else { + if (RES_OK == setup_next_segment(&rs)) { propagate(&rs); } - } while (rs.end == TERM_NONE); + } while (!rs.end); } - /* propagation ended */ -#if 0 - /* TODO: use this to feed the implicit MC computations */ - fprintf(output, "Realisation %u success mask: 0x%0x\n", - (unsigned) r, rs.success_mask); - fprintf(output, "Realisation %u end: %s\n\n", - (unsigned) r, END_TEXT[rs.end]); -#endif // 0 - - /* must retry failed realisations */ - if (rs.end == TERM_ERR) { + if (rs.error) { + /* TODO: must retry failed realisations */ log_error(scene->dev, "%s: realisation failure.\n", FUNC_NAME); goto error; } - else r++; - } while (rs.end == TERM_ERR); + /* propagation ended: feed implicit MC data */ + seg = current_segment(&rs); + w = seg->weight; + + ASSERT(!rs.success | !rs.shadow); + if (rs.shadow) { + estimator->shadow.weight += w; + estimator->shadow.sqr_weight += w * w; + } + else if (!rs.success) { + estimator->missing.weight += w; + estimator->missing.sqr_weight += w * w; + } + + r++; + } while (rs.error); } + estimator->realisation_count += realisations_count; exit: release_realisation(&rs); return res; diff --git a/src/ssol_solver_c.h b/src/ssol_solver_c.h @@ -42,26 +42,6 @@ struct ranst_sun_wl; #define DARRAY_NAME instances_ptr #include <rsys/dynamic_array.h> -enum realisation_termination { - TERM_NONE, - TERM_SUCCESS, - TERM_SHADOW, - TERM_POINTING, - TERM_MISSING, - TERM_BLOCKED, - TERM_ERR, - - TERM_COUNT__ -}; - -enum realisation_mode { - MODE_NONE, - MODE_STD, - MODE_ROULETTE, - - MODE_COUNT__ -}; - struct segment { const struct ssol_instance* hit_instance; const struct ssol_instance* self_instance; /* instance of the starting point */ @@ -132,15 +112,14 @@ struct solver_data { }; struct realisation { - enum realisation_termination end; - enum realisation_mode mode; struct darray_segment segments; struct starting_point start; struct solver_data data; double wavelength; size_t rs_id; uint32_t s_idx; - uint32_t success_mask; + /* status */ + char end, error, shadow, success; }; extern LOCAL_SYM res_T diff --git a/src/test_ssol_solver1.c b/src/test_ssol_solver1.c @@ -56,6 +56,8 @@ main(int argc, char** argv) struct ssol_spectrum* spectrum; struct ssol_spectrum* abs; struct ssol_atmosphere* atm; + struct ssol_estimator* estimator; + struct ssol_estimator_status status; double dir[3]; double wavelengths[3] = { 1, 2, 3 }; double mismatch[3] = { 1.5, 3.5 }; @@ -100,12 +102,14 @@ main(int argc, char** argv) CHECK(ssol_sun_set_dni(sun, 1000), RES_OK); CHECK(ssol_scene_create(dev, &scene), RES_OK); CHECK(ssol_scene_attach_sun(scene, sun), RES_OK); + CHECK(ssol_estimator_create(dev, &estimator), RES_OK); - CHECK(ssol_solve(NULL, rng, 10, stdout), RES_BAD_ARG); - CHECK(ssol_solve(scene, NULL, 10, stdout), RES_BAD_ARG); - CHECK(ssol_solve(scene, rng, 0, stdout), RES_BAD_ARG); - CHECK(ssol_solve(scene, rng, 10, NULL), RES_BAD_ARG); - CHECK(ssol_solve(scene, rng, 10, stdout), RES_BAD_ARG); /* no geometry */ + CHECK(ssol_solve(NULL, rng, 10, stdout, estimator), RES_BAD_ARG); + CHECK(ssol_solve(scene, NULL, 10, stdout, estimator), RES_BAD_ARG); + CHECK(ssol_solve(scene, rng, 0, stdout, estimator), RES_BAD_ARG); + CHECK(ssol_solve(scene, rng, 10, NULL, estimator), RES_BAD_ARG); + CHECK(ssol_solve(scene, rng, 10, stdout, NULL), RES_BAD_ARG); + CHECK(ssol_solve(scene, rng, 10, stdout, estimator), RES_BAD_ARG); /* no geometry */ /* create scene content */ @@ -141,25 +145,25 @@ main(int argc, char** argv) CHECK(ssol_instance_set_target_mask(target, 0x4, 0), RES_OK); CHECK(ssol_scene_attach_instance(scene, target), RES_OK); - CHECK(ssol_solve(scene, rng, 1, stdout), RES_OK); /* ready to solve! */ + CHECK(ssol_solve(scene, rng, 1, stdout, estimator), RES_OK); /* ready to solve! */ CHECK(ssol_instance_dont_sample(target, 1), RES_OK); CHECK(ssol_instance_dont_sample(secondary, 1), RES_OK); CHECK(ssol_instance_dont_sample(heliostat, 1), RES_OK); - CHECK(ssol_solve(scene, rng, 10, stdout), RES_BAD_ARG); /* no geometry to sample */ + CHECK(ssol_solve(scene, rng, 10, stdout, estimator), RES_BAD_ARG); /* no geometry to sample */ CHECK(ssol_instance_dont_sample(target, 0), RES_OK); CHECK(ssol_instance_dont_sample(secondary, 0), RES_OK); CHECK(ssol_instance_dont_sample(heliostat, 0), RES_OK); CHECK(ssol_scene_detach_sun(scene, sun), RES_OK); - CHECK(ssol_solve(scene, rng, 10, stdout), RES_BAD_ARG); /* no attached sun */ + CHECK(ssol_solve(scene, rng, 10, stdout, estimator), RES_BAD_ARG); /* no attached sun */ CHECK(ssol_sun_ref_put(sun), RES_OK); CHECK(ssol_sun_create_directional(dev, &sun), RES_OK); 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, stdout), RES_BAD_ARG); /* sun with no spectrum */ + CHECK(ssol_solve(scene, rng, 10, stdout, estimator), RES_BAD_ARG); /* sun with no spectrum */ CHECK(ssol_scene_detach_sun(scene, sun), RES_OK); CHECK(ssol_sun_ref_put(sun), RES_OK); @@ -167,13 +171,13 @@ 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, stdout), RES_BAD_ARG); /* sun with undefined DNI */ + CHECK(ssol_solve(scene, rng, 10, stdout, estimator), RES_BAD_ARG); /* sun with undefined DNI */ CHECK(ssol_sun_set_dni(sun, 1000), RES_OK); CHECK(ssol_instance_set_receiver(heliostat, NULL, NULL), RES_OK); CHECK(ssol_instance_set_receiver(secondary, NULL, NULL), RES_OK); CHECK(ssol_instance_set_receiver(target, NULL, NULL), RES_OK); - CHECK(ssol_solve(scene, rng, 10, stdout), RES_BAD_ARG); /* no receiver in scene */ + CHECK(ssol_solve(scene, rng, 10, stdout, estimator), RES_BAD_ARG); /* no receiver in scene */ CHECK(ssol_instance_set_receiver(heliostat, "miroir", NULL), RES_OK); CHECK(ssol_instance_set_receiver(secondary, "secondaire", NULL), RES_OK); CHECK(ssol_instance_set_receiver(target, "cible", NULL), RES_OK); @@ -183,25 +187,33 @@ main(int argc, char** argv) CHECK(ssol_atmosphere_create_uniform(dev, &atm), RES_OK); CHECK(ssol_atmosphere_set_uniform_absorbtion(atm, abs), RES_OK); CHECK(ssol_scene_attach_atmosphere(scene, atm), RES_OK); - CHECK(ssol_solve(scene, rng, 10, stdout), RES_BAD_ARG); /* spectra mismatch */ + CHECK(ssol_solve(scene, rng, 10, stdout, estimator), RES_BAD_ARG); /* spectra mismatch */ 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); /* can sample any geometry; variance is high */ - tmp = tmpfile(); - CHECK(!tmp, 0); + NCHECK(tmp = tmpfile(), 0); #define N 5000 - CHECK(ssol_solve(scene, rng, N, tmp), RES_OK); + CHECK(ssol_estimator_clear(estimator), RES_OK); + CHECK(ssol_solve(scene, rng, N, tmp, estimator), RES_OK); CHECK(get_receiver_id(target, 1, &r_id), RES_OK); CHECK(pp_sum(tmp, r_id, N, &m, &std), RES_OK); - logger_print(&logger, LOG_OUTPUT, "\nP = %g +/- %g\n", m, std); -#define DNI_cos (1000 * cos(PI / 4)) + logger_print(&logger, LOG_OUTPUT, "\nP = %g +/- %g", m, std); +#define COS cos(PI / 4) +#define DNI_cos (1000 * COS) CHECK(eq_eps(m, 4 * DNI_cos, MMAX(4 * DNI_cos * 1e-2, std)), 1); #define SQR(x) ((x)*(x)) dbl = sqrt((SQR(12 * DNI_cos) / 3 - SQR(4 * DNI_cos)) / N); CHECK(eq_eps(std, dbl, dbl*1e-2), 1); CHECK(fclose(tmp), 0); + /* target was sampled but shadowed by secondary */ + CHECK(ssol_estimator_get_status(estimator, STATUS_SHADOW, &status), RES_OK); + logger_print(&logger, LOG_OUTPUT, "Shadows = %g +/- %g", status.E, status.SE); + CHECK(eq_eps(status.E, m, 2 * dbl), 1); + CHECK(ssol_estimator_get_status(estimator, STATUS_MISSING, &status), RES_OK); + logger_print(&logger, LOG_OUTPUT, "Missing = %g +/- %g", status.E, status.SE); + CHECK(eq_eps(status.E, m, status.SE), 1); /* sample primary mirror only; variance is low */ CHECK(ssol_instance_dont_sample(target, 1), RES_OK); @@ -210,13 +222,20 @@ main(int argc, char** argv) CHECK(ssol_instance_set_target_mask(secondary, 0, 0), RES_OK); CHECK(ssol_instance_set_target_mask(target, 0x1, 0), RES_OK); - tmp = tmpfile(); - CHECK(ssol_solve(scene, rng, N, tmp), RES_OK); + NCHECK(tmp = tmpfile(), 0); + CHECK(ssol_estimator_clear(estimator), RES_OK); + CHECK(ssol_solve(scene, rng, N, tmp, estimator), RES_OK); CHECK(pp_sum(tmp, r_id, N, &m, &std), RES_OK); - logger_print(&logger, LOG_OUTPUT, "\nP = %g +/- %g\n", m, std); + logger_print(&logger, LOG_OUTPUT, "\nP = %g +/- %g", m, std); CHECK(eq_eps(m, 4 * DNI_cos, MMAX(4 * DNI_cos * 1e-2, std)), 1); CHECK(eq_eps(std, 0, 1e-4), 1); CHECK(fclose(tmp), 0); + CHECK(ssol_estimator_get_status(estimator, STATUS_SHADOW, &status), RES_OK); + logger_print(&logger, LOG_OUTPUT, "Shadows = %g +/- %g", status.E, status.SE); + CHECK(eq_eps(status.E, 0, 1e-4), 1); + CHECK(ssol_estimator_get_status(estimator, STATUS_MISSING, &status), RES_OK); + logger_print(&logger, LOG_OUTPUT, "Missing = %g +/- %g", status.E, status.SE); + CHECK(eq_eps(status.E, 0, 1e-4), 1); /* check atmosphere model; with no absorbtion result is unchanged */ CHECK(ssol_spectrum_create(dev, &abs), RES_OK); @@ -226,15 +245,22 @@ 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), RES_OK); + CHECK(ssol_estimator_clear(estimator), RES_OK); + CHECK(ssol_solve(scene, rng, N, tmp, estimator), RES_OK); CHECK(pp_sum(tmp, r_id, N, &m, &std), RES_OK); - logger_print(&logger, LOG_OUTPUT, "\nP = %g +/- %g\n", m, std); + logger_print(&logger, LOG_OUTPUT, "\nP = %g +/- %g", m, std); CHECK(fclose(tmp), 0); CHECK(eq_eps(m, 4 * DNI_cos, MMAX(4 * DNI_cos * 1e-2, std)), 1); CHECK(eq_eps(std, 0, 1e-4), 1); 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); + CHECK(ssol_estimator_get_status(estimator, STATUS_SHADOW, &status), RES_OK); + logger_print(&logger, LOG_OUTPUT, "Shadows = %g +/- %g", status.E, status.SE); + CHECK(eq_eps(status.E, 0, 1e-4), 1); + CHECK(ssol_estimator_get_status(estimator, STATUS_MISSING, &status), RES_OK); + logger_print(&logger, LOG_OUTPUT, "Missing = %g +/- %g", status.E, status.SE); + CHECK(eq_eps(status.E, 0, 1e-4), 1); /* check atmosphere model; with absorbtion power decreases */ ka[0] = ka[1] = ka[2] = 0.1; @@ -245,13 +271,20 @@ 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), RES_OK); + CHECK(ssol_estimator_clear(estimator), RES_OK); + CHECK(ssol_solve(scene, rng, N, tmp, estimator), RES_OK); CHECK(pp_sum(tmp, r_id, N, &m, &std), RES_OK); - logger_print(&logger, LOG_OUTPUT, "\nP = %g +/- %g\n", m, std); + logger_print(&logger, LOG_OUTPUT, "\nP = %g +/- %g", m, std); CHECK(fclose(tmp), 0); #define K (exp(-0.1 * 4 * sqrt(2))) CHECK(eq_eps(m, 4 * K * DNI_cos, MMAX(4 * K * DNI_cos * 1e-1, std)), 1); CHECK(eq_eps(std, 0, 1e-4), 1); + CHECK(ssol_estimator_get_status(estimator, STATUS_SHADOW, &status), RES_OK); + logger_print(&logger, LOG_OUTPUT, "Shadows = %g +/- %g", status.E, status.SE); + CHECK(eq_eps(status.E, 0, 1e-4), 1); + CHECK(ssol_estimator_get_status(estimator, STATUS_MISSING, &status), RES_OK); + logger_print(&logger, LOG_OUTPUT, "Missing = %g +/- %g", status.E, status.SE); + CHECK(eq_eps(status.E, 0, 1e-4), 1); /* check a monochromatic sun */ CHECK(ssol_spectrum_setup(spectrum, &mono, intensities, 1), RES_OK); @@ -264,13 +297,20 @@ main(int argc, char** argv) ka[1] = 0.2; CHECK(ssol_spectrum_setup(abs, wavelengths, ka, 2), RES_OK); NCHECK(tmp = tmpfile(), 0); - CHECK(ssol_solve(scene, rng, N, tmp), RES_OK); + CHECK(ssol_estimator_clear(estimator), RES_OK); + CHECK(ssol_solve(scene, rng, N, tmp, estimator), RES_OK); CHECK(pp_sum(tmp, r_id, N, &m, &std), RES_OK); - logger_print(&logger, LOG_OUTPUT, "\nP = %g +/- %g\n", m, std); + logger_print(&logger, LOG_OUTPUT, "\nP = %g +/- %g", m, std); CHECK(fclose(tmp), 0); #define K2 (exp(-0.121 * 4 * sqrt(2))) CHECK(eq_eps(m, 4 * K2 * DNI_cos, MMAX(4 * K2 * DNI_cos * 1e-4, std)), 1); CHECK(eq_eps(std, 0, 1e-4), 1); + CHECK(ssol_estimator_get_status(estimator, STATUS_SHADOW, &status), RES_OK); + logger_print(&logger, LOG_OUTPUT, "Shadows = %g +/- %g", status.E, status.SE); + CHECK(eq_eps(status.E, 0, 1e-4), 1); + CHECK(ssol_estimator_get_status(estimator, STATUS_MISSING, &status), RES_OK); + logger_print(&logger, LOG_OUTPUT, "Missing = %g +/- %g", status.E, status.SE); + CHECK(eq_eps(status.E, 0, 1e-4), 1); /* free data */ @@ -287,6 +327,7 @@ main(int argc, char** argv) CHECK(ssp_rng_ref_put(rng), RES_OK); CHECK(ssol_spectrum_ref_put(abs), RES_OK); CHECK(ssol_atmosphere_ref_put(atm), RES_OK); + CHECK(ssol_estimator_ref_put(estimator), RES_OK); CHECK(ssol_spectrum_ref_put(spectrum), RES_OK); CHECK(ssol_sun_ref_put(sun), RES_OK); CHECK(ssol_sun_ref_put(sun_mono), RES_OK); diff --git a/src/test_ssol_solver2.c b/src/test_ssol_solver2.c @@ -70,6 +70,8 @@ main(int argc, char** argv) struct ssol_instance* target; struct ssol_sun* sun; struct ssol_spectrum* spectrum; + struct ssol_estimator* estimator; + struct ssol_estimator_status status; double dir[3]; double wavelengths[3] = { 1, 2, 3 }; double intensities[3] = { 1, 0.8, 1 }; @@ -114,6 +116,7 @@ main(int argc, char** argv) CHECK(ssol_sun_set_dni(sun, 1000), RES_OK); CHECK(ssol_scene_create(dev, &scene), RES_OK); CHECK(ssol_scene_attach_sun(scene, sun), RES_OK); + CHECK(ssol_estimator_create(dev, &estimator), RES_OK); /* create scene content */ @@ -174,18 +177,24 @@ main(int argc, char** argv) CHECK(ssol_instance_dont_sample(target, 1), RES_OK); CHECK(ssol_scene_attach_instance(scene, target), RES_OK); - CHECK(ssol_solve(scene, rng, 20, stdout), RES_OK); - - tmp = tmpfile(); + NCHECK(tmp = tmpfile(), 0); #define N 10000 - CHECK(ssol_solve(scene, rng, N, tmp), RES_OK); + CHECK(ssol_solve(scene, rng, N, tmp, estimator), RES_OK); CHECK(get_receiver_id(target, 1, &r_id), RES_OK); CHECK(pp_sum(tmp, r_id, N, &m, &std), RES_OK); + CHECK(fclose(tmp), 0); logger_print(&logger, LOG_OUTPUT, "\nP = %g +/- %g\n", m, std); #define DNI_cos (1000 * cos(PI / 4)) 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_status(estimator, STATUS_SHADOW, &status), RES_OK); + logger_print(&logger, LOG_OUTPUT, "Shadows = %g +/- %g", status.E, status.SE); + CHECK(eq_eps(status.E, 0, 1e-4), 1); + CHECK(ssol_estimator_get_status(estimator, STATUS_MISSING, &status), RES_OK); + logger_print(&logger, LOG_OUTPUT, "Missing = %g +/- %g", status.E, status.SE); + CHECK(eq_eps(status.E, 0, 1e-4), 1); + /* free data */ CHECK(ssol_instance_ref_put(heliostat1), RES_OK); @@ -201,6 +210,7 @@ main(int argc, char** argv) CHECK(ssol_material_ref_put(m_mtl), RES_OK); CHECK(ssol_material_ref_put(v_mtl), RES_OK); CHECK(ssol_device_ref_put(dev), RES_OK); + CHECK(ssol_estimator_ref_put(estimator), 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); diff --git a/src/test_ssol_solver2b.c b/src/test_ssol_solver2b.c @@ -70,6 +70,8 @@ main(int argc, char** argv) struct ssol_instance* target; struct ssol_sun* sun; struct ssol_spectrum* spectrum; + struct ssol_estimator* estimator; + struct ssol_estimator_status status; double dir[3]; double wavelengths[3] = { 1, 2, 3 }; double intensities[3] = { 1, 0.8, 1 }; @@ -114,6 +116,7 @@ main(int argc, char** argv) CHECK(ssol_sun_set_dni(sun, 1000), RES_OK); CHECK(ssol_scene_create(dev, &scene), RES_OK); CHECK(ssol_scene_attach_sun(scene, sun), RES_OK); + CHECK(ssol_estimator_create(dev, &estimator), RES_OK); /* create scene content */ @@ -178,19 +181,25 @@ main(int argc, char** argv) CHECK(ssol_instance_set_target_mask(target, 0x1, 0), RES_OK); CHECK(ssol_instance_dont_sample(target, 1), RES_OK); CHECK(ssol_scene_attach_instance(scene, target), RES_OK); - - CHECK(ssol_solve(scene, rng, 20, stdout), RES_OK); - - tmp = tmpfile(); + + NCHECK(tmp = tmpfile(), 0); #define N 50000 - CHECK(ssol_solve(scene, rng, N, tmp), RES_OK); + CHECK(ssol_solve(scene, rng, N, tmp, estimator), RES_OK); CHECK(get_receiver_id(target, 1, &r_id), RES_OK); CHECK(pp_sum(tmp, r_id, N, &m, &std), RES_OK); + CHECK(fclose(tmp), 0); logger_print(&logger, LOG_OUTPUT, "\nP = %g +/- %g\n", m, std); #define DNI_cos (1000 * cos(PI / 4)) CHECK(eq_eps(m, 2 * DNI_cos, MMAX(2 * DNI_cos * 1e-2, std)), 1); #define SQR(x) ((x)*(x)) CHECK(eq_eps(std, sqrt((SQR(4 * DNI_cos) / 2 - SQR(2 * DNI_cos)) / N), 1e-3), 1); + CHECK(ssol_estimator_get_status(estimator, STATUS_SHADOW, &status), RES_OK); + logger_print(&logger, LOG_OUTPUT, "Shadows = %g +/- %g", status.E, status.SE); + CHECK(eq_eps(status.E, 0, 1e-4), 1); + CHECK(ssol_estimator_get_status(estimator, STATUS_MISSING, &status), RES_OK); + logger_print(&logger, LOG_OUTPUT, "Missing = %g +/- %g", status.E, status.SE); + CHECK(eq_eps(status.E, 0, 1e-4), 1); + /* free data */ CHECK(ssol_instance_ref_put(heliostat1), RES_OK); @@ -205,6 +214,7 @@ main(int argc, char** argv) 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); diff --git a/src/test_ssol_solver3.c b/src/test_ssol_solver3.c @@ -61,6 +61,8 @@ main(int argc, char** argv) struct ssol_instance* target; struct ssol_sun* sun; struct ssol_spectrum* spectrum; + struct ssol_estimator* estimator; + struct ssol_estimator_status status; double dir[3]; double wavelengths[3] = { 1, 2, 3 }; double intensities[3] = { 1, 0.8, 1 }; @@ -94,6 +96,7 @@ main(int argc, char** argv) CHECK(ssol_sun_set_dni(sun, 1000), RES_OK); CHECK(ssol_scene_create(dev, &scene), RES_OK); CHECK(ssol_scene_attach_sun(scene, sun), RES_OK); + CHECK(ssol_estimator_create(dev, &estimator), RES_OK); /* create scene content */ @@ -133,17 +136,24 @@ main(int argc, char** argv) CHECK(ssol_instance_set_target_mask(target, 0x1, 0), RES_OK); CHECK(ssol_instance_dont_sample(target, 1), RES_OK); CHECK(ssol_scene_attach_instance(scene, target), RES_OK); - - tmp = tmpfile(); + + NCHECK(tmp = tmpfile(), 0); #define N 20000 - CHECK(ssol_solve(scene, rng, N, tmp), RES_OK); + CHECK(ssol_solve(scene, rng, N, tmp, estimator), RES_OK); CHECK(get_receiver_id(target, 1, &r_id), RES_OK); CHECK(pp_sum(tmp, r_id, N, &m, &std), RES_OK); + CHECK(fclose(tmp), 0); logger_print(&logger, LOG_OUTPUT, "\nP = %g +/- %g\n", m, std); #define DNI_cos (1000 * cos(PI / 4)) CHECK(eq_eps(m, 4 * DNI_cos, 4 * DNI_cos * 2e-1), 1); #define SQR(x) ((x)*(x)) CHECK(eq_eps(std, sqrt((SQR(400 * DNI_cos) / 100 - SQR(4 * DNI_cos)) / N), 20), 1); + CHECK(ssol_estimator_get_status(estimator, STATUS_SHADOW, &status), RES_OK); + logger_print(&logger, LOG_OUTPUT, "Shadows = %g +/- %g", status.E, status.SE); + CHECK(eq_eps(status.E, 0, 1e-4), 1); + CHECK(ssol_estimator_get_status(estimator, STATUS_MISSING, &status), RES_OK); + logger_print(&logger, LOG_OUTPUT, "Missing = %g +/- %g", status.E, status.SE); + CHECK(eq_eps(status.E, 0, 1e-4), 1); /* free data */ @@ -156,6 +166,7 @@ main(int argc, char** argv) CHECK(ssol_material_ref_put(m_mtl), RES_OK); CHECK(ssol_material_ref_put(v_mtl), RES_OK); CHECK(ssol_device_ref_put(dev), RES_OK); + CHECK(ssol_estimator_ref_put(estimator), 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); diff --git a/src/test_ssol_solver3N.c b/src/test_ssol_solver3N.c @@ -104,6 +104,8 @@ main(int argc, char** argv) struct ssol_instance* target; struct ssol_sun* sun; struct ssol_spectrum* spectrum; + struct ssol_estimator* estimator; + struct ssol_estimator_status status; double sun_dir[3]; double target_pos[3]; double wavelengths[3] = { 1, 2, 3 }; @@ -136,6 +138,7 @@ main(int argc, char** argv) CHECK(ssol_sun_set_dni(sun, 1000), RES_OK); CHECK(ssol_scene_create(dev, &scene), RES_OK); CHECK(ssol_scene_attach_sun(scene, sun), RES_OK); + CHECK(ssol_estimator_create(dev, &estimator), RES_OK); /* create scene content */ @@ -199,14 +202,21 @@ main(int argc, char** argv) CHECK(ssol_instance_dont_sample(target, 1), RES_OK); CHECK(ssol_scene_attach_instance(scene, target), RES_OK); - tmp = tmpfile(); + NCHECK(tmp = tmpfile(), 0); #define N 10000 - CHECK(ssol_solve(scene, rng, N, tmp), RES_OK); + CHECK(ssol_solve(scene, rng, N, tmp, estimator), RES_OK); CHECK(get_receiver_id(target, 1, &r_id), RES_OK); CHECK(pp_sum(tmp, r_id, N, &m, &std), RES_OK); + CHECK(fclose(tmp), 0); logger_print(&logger, LOG_OUTPUT, "\nP = %g +/- %g\n", m, std); #define DNI_cos (1000 * cos(PI / 8)) CHECK(eq_eps(m, 4 * NX * NY * NZ * DNI_cos, 4 * NX * NY * NZ * DNI_cos * 2e-1), 1); + CHECK(ssol_estimator_get_status(estimator, STATUS_SHADOW, &status), RES_OK); + logger_print(&logger, LOG_OUTPUT, "Shadows = %g +/- %g", status.E, status.SE); + CHECK(eq_eps(status.E, 0, 1e-4), 1); + CHECK(ssol_estimator_get_status(estimator, STATUS_MISSING, &status), RES_OK); + logger_print(&logger, LOG_OUTPUT, "Missing = %g +/- %g", status.E, status.SE); + CHECK(eq_eps(status.E, 0, 1e-4), 1); /* free data */ @@ -217,6 +227,7 @@ main(int argc, char** argv) 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); diff --git a/src/test_ssol_solver4.c b/src/test_ssol_solver4.c @@ -62,6 +62,8 @@ main(int argc, char** argv) struct ssol_instance* target2; struct ssol_sun* sun; struct ssol_spectrum* spectrum; + struct ssol_estimator* estimator; + struct ssol_estimator_status status; double dir[3]; double wavelengths[3] = { 1, 2, 3 }; double intensities[3] = { 1, 0.8, 1 }; @@ -95,6 +97,7 @@ main(int argc, char** argv) CHECK(ssol_sun_set_dni(sun, 1000), RES_OK); CHECK(ssol_scene_create(dev, &scene), RES_OK); CHECK(ssol_scene_attach_sun(scene, sun), RES_OK); + CHECK(ssol_estimator_create(dev, &estimator), RES_OK); /* create scene content */ @@ -140,20 +143,27 @@ main(int argc, char** argv) CHECK(ssol_instance_set_target_mask(target2, 0x1, 0), RES_OK); CHECK(ssol_instance_dont_sample(target2, 1), RES_OK); CHECK(ssol_scene_attach_instance(scene, target2), RES_OK); - - tmp = tmpfile(); + + NCHECK(tmp = tmpfile(), 0); #define N 10000 - CHECK(ssol_solve(scene, rng, N, tmp), RES_OK); + CHECK(ssol_solve(scene, rng, N, tmp, estimator), RES_OK); CHECK(get_receiver_id(target1, 1, &r_id1), RES_OK); CHECK(get_receiver_id(target2, 1, &r_id2), RES_OK); CHECK(pp_sum(tmp, r_id1, N, &m1, &std1), RES_OK); CHECK(pp_sum(tmp, r_id2, N, &m2, &std2), RES_OK); + CHECK(fclose(tmp), 0); logger_print(&logger, LOG_OUTPUT, "\nP = %g +/- %g\n", m1, std1); #define DNI_cos (1000 * cos(0)) CHECK(eq_eps(m1, 400 * DNI_cos, 400 * DNI_cos * 1e-4), 1); CHECK(eq_eps(std1, 0, 1), 1); CHECK(m1, m2); CHECK(std1, std2); + CHECK(ssol_estimator_get_status(estimator, STATUS_SHADOW, &status), RES_OK); + logger_print(&logger, LOG_OUTPUT, "Shadows = %g +/- %g", status.E, status.SE); + CHECK(eq_eps(status.E, 0, 1e-4), 1); + CHECK(ssol_estimator_get_status(estimator, STATUS_MISSING, &status), RES_OK); + logger_print(&logger, LOG_OUTPUT, "Missing = %g +/- %g", status.E, status.SE); + CHECK(eq_eps(status.E, 0, 1e-4), 1); /* free data */ @@ -166,6 +176,7 @@ main(int argc, char** argv) 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); diff --git a/src/test_ssol_solver5.c b/src/test_ssol_solver5.c @@ -61,6 +61,8 @@ main(int argc, char** argv) struct ssol_instance* target; struct ssol_sun* sun; struct ssol_spectrum* spectrum; + struct ssol_estimator* estimator; + struct ssol_estimator_status status; double dir[3]; double wavelengths[3] = { 1, 2, 3 }; double intensities[3] = { 1, 0.8, 1 }; @@ -94,6 +96,7 @@ main(int argc, char** argv) CHECK(ssol_sun_set_dni(sun, 1000), RES_OK); CHECK(ssol_scene_create(dev, &scene), RES_OK); CHECK(ssol_scene_attach_sun(scene, sun), RES_OK); + CHECK(ssol_estimator_create(dev, &estimator), RES_OK); /* create scene content */ @@ -134,15 +137,22 @@ main(int argc, char** argv) CHECK(ssol_instance_dont_sample(target, 1), RES_OK); CHECK(ssol_scene_attach_instance(scene, target), RES_OK); - tmp = tmpfile(); + NCHECK(tmp = tmpfile(), 0); #define N 10000 - CHECK(ssol_solve(scene, rng, N, tmp), RES_OK); + CHECK(ssol_solve(scene, rng, N, tmp, estimator), RES_OK); CHECK(get_receiver_id(target, 1, &r_id), RES_OK); CHECK(pp_sum(tmp, r_id, N, &m, &std), RES_OK); + CHECK(fclose(tmp), 0); logger_print(&logger, LOG_OUTPUT, "\nP = %g +/- %g\n", m, std); #define DNI_cos (1000 * cos(0)) CHECK(eq_eps(m, 400 * DNI_cos, 20), 1); CHECK(eq_eps(std, 0, 1), 1); + CHECK(ssol_estimator_get_status(estimator, STATUS_SHADOW, &status), RES_OK); + logger_print(&logger, LOG_OUTPUT, "Shadows = %g +/- %g", status.E, status.SE); + CHECK(eq_eps(status.E, 0, 1e-4), 1); + CHECK(ssol_estimator_get_status(estimator, STATUS_MISSING, &status), RES_OK); + logger_print(&logger, LOG_OUTPUT, "Missing = %g +/- %g", status.E, status.SE); + CHECK(eq_eps(status.E, 0, 1e-4), 1); /* free data */ @@ -152,6 +162,7 @@ main(int argc, char** argv) CHECK(ssol_object_ref_put(t_object), RES_OK); CHECK(ssol_shape_ref_put(rect), RES_OK); CHECK(ssol_shape_ref_put(quad_square), RES_OK); + CHECK(ssol_estimator_ref_put(estimator), RES_OK); CHECK(ssol_material_ref_put(m_mtl), RES_OK); CHECK(ssol_material_ref_put(v_mtl), RES_OK); CHECK(ssol_device_ref_put(dev), RES_OK);