solstice-solver

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

commit 5759956711e0c67044f749ff2bc762cf3b89d5a3
parent 937d670760f5ce6c0829c2a4f47b266a76335e13
Author: Christophe Coustet <christophe.coustet@meso-star.com>
Date:   Tue,  6 Sep 2016 16:47:58 +0200

Add some code for quadrics.

Recompute hits on quadrics; both position and normal.
Need to add world/local transforms.

Diffstat:
Msrc/ssol_instance_c.h | 7+++++++
Msrc/ssol_scene.c | 56++++++++++++++++++++++++++++++++++++++++++--------------
Msrc/ssol_shape.c | 145+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Msrc/ssol_shape_c.h | 41+++++++++++++++++++++++++++++++++++++++++
Msrc/ssol_solver.c | 32+++-----------------------------
Msrc/ssol_solver_c.h | 1+
Msrc/ssol_sun_c.h | 7-------
7 files changed, 239 insertions(+), 50 deletions(-)

diff --git a/src/ssol_instance_c.h b/src/ssol_instance_c.h @@ -41,4 +41,11 @@ instance_get_receiver_name(const struct ssol_instance* instance) ? NULL : str_cget(&instance->receiver_name); } +static INLINE int +is_instance_punched(const struct ssol_instance* instance) +{ + ASSERT(instance); + return instance->object->shape->type == SHAPE_PUNCHED; +} + #endif /* SSOL_INSTANCE_C_H */ diff --git a/src/ssol_scene.c b/src/ssol_scene.c @@ -20,13 +20,14 @@ #include "ssol_sun_c.h" #include "ssol_device_c.h" #include "ssol_material_c.h" +#include "ssol_shape_c.h" #include "ssol_object_c.h" #include "ssol_instance_c.h" #include <rsys/list.h> #include <rsys/mem_allocator.h> #include <rsys/rsys.h> -#include <rsys/float3.h> +#include <rsys/double3.h> /******************************************************************************* * Helper functions @@ -115,7 +116,7 @@ ssol_scene_attach_instance res_T res; if(!scene || !instance) return RES_BAD_ARG; - + /* Attach the instantiated s3d shape to ray-trace to the RT scene */ res = s3d_scene_attach_shape(scene->scn_rt, instance->shape_rt); if(res != RES_OK) return res; @@ -229,24 +230,24 @@ scene_setup_s3d_sampling_scene(struct ssol_scene* scn) htable_instance_begin(&scn->instances_rt, &it); htable_instance_end(&scn->instances_rt, &end); - while(!htable_instance_iterator_eq(&it, &end)) { + while (!htable_instance_iterator_eq(&it, &end)) { struct ssol_instance* inst = *htable_instance_iterator_data_get(&it); unsigned id; htable_instance_iterator_next(&it); /* TODO: keep only primary mirrors */ - if(inst->object->material->type != MATERIAL_MIRROR) + if (inst->object->material->type != MATERIAL_MIRROR) continue; /* Attach the instantiated s3d sampling shape to the s3d sampling scene */ res = s3d_scene_attach_shape(scn->scn_samp, inst->shape_samp); - if(res != RES_OK) goto error; + if (res != RES_OK) goto error; /* Register the instantiated s3d sampling shape */ S3D(shape_get_id(inst->shape_samp, &id)); ASSERT(!htable_instance_find(&scn->instances_samp, &id)); res = htable_instance_set(&scn->instances_samp, &id, &inst); - if(res != RES_OK) goto error; + if (res != RES_OK) goto error; /* Do not get a reference onto the instance since it was already referenced * by the scene on its attachment */ @@ -278,7 +279,7 @@ hit_filter_function struct segment* prev; int front_face = 0; - (void) filter_data; + (void) filter_data, (void) org, (void) dir; ASSERT(rs); prev = previous_segment(rs); seg = current_segment(rs); @@ -290,17 +291,44 @@ hit_filter_function return 1; /* Discard self intersection */ inst = *htable_instance_find(&rs->data.scene->instances_rt, &hit->prim.inst_id); + + if (is_instance_punched(inst)) { + /* hits on quadrics must be recomputed more accurately */ + double dist; + switch (inst->object->shape->quadric.type) { + case SSOL_QUADRIC_PLANE: { + int ok = quadric_plane_intersect_local( + seg->org, seg->dir, seg->hit_pos, seg->hit_normal, &dist); + if (!ok) return 1; /* invalid impact */ + break; + } + case SSOL_QUADRIC_PARABOLIC_CYLINDER: { + const struct ssol_quadric_parabolic_cylinder* quad + = (struct ssol_quadric_parabolic_cylinder*)&inst->object->shape->quadric; + int ok = quadric_parabolic_cylinder_intersect_local( + seg->org, seg->dir, quad, seg->hit_pos, seg->hit_normal, &dist); + if (!ok) return 1; /* invalid impact */ + break; + } + case SSOL_QUADRIC_PARABOL: { + const struct ssol_quadric_parabol* quad + = (struct ssol_quadric_parabol*)&inst->object->shape->quadric; + int ok = quadric_parabol_intersect_local( + seg->org, seg->dir, quad, seg->hit_pos, seg->hit_normal, &dist); + if (!ok) return 1; /* invalid impact */ + break; + } + default: FATAL("Unreachable code\n"); break; + } + } /* Check if the hit surface is a receiver that registers hit data */ receiver_name = instance_get_receiver_name(inst); if(receiver_name) { - float cos_in; + double cos_in; /* check normal orientation */ - cos_in = f3_dot(hit->normal, dir); + cos_in = d3_dot(seg->hit_normal, seg->dir); if (cos_in < 0) { - float pos[3]; - - f3_set(pos, f3_add(pos, org, f3_mulf(pos, dir, hit->distance))); front_face = 1; fprintf(rs->data.out_stream, "Receiver '%s': %u %u %g %g (%g:%g:%g) (%g:%g:%g) (%g:%g)\n", @@ -309,8 +337,8 @@ hit_filter_function (unsigned)rs->s_idx, rs->freq, seg->weight, - SPLIT3(pos), - SPLIT3(dir), + SPLIT3(seg->hit_pos), + SPLIT3(seg->dir), SPLIT2(hit->uv)); } } diff --git a/src/ssol_shape.c b/src/ssol_shape.c @@ -18,6 +18,7 @@ #include "ssol_device_c.h" #include "ssol_shape_c.h" +#include <rsys/double3.h> #include <rsys/double2.h> #include <rsys/dynamic_array_double.h> #include <rsys/dynamic_array_size_t.h> @@ -537,6 +538,150 @@ shape_release(ref_T* ref) } /******************************************************************************* +* Local functions +******************************************************************************/ +void +quadric_plane_gradient_local(double grad[3]) +{ + grad[0] = 0; + grad[1] = 0; + grad[2] = 1; +} + +void +quadric_parabol_gradient_local + (double grad[3], + const double pt[3], + const struct ssol_quadric_parabol* quad) +{ + grad[0] = -pt[0]; + grad[1] = -pt[1]; + grad[2] = 2 * quad->focal; +} + +void +quadric_parabolic_cylinder_gradient_local + (double grad[3], + const double pt[3], + const struct ssol_quadric_parabolic_cylinder* quad) +{ + grad[0] = 0; + grad[1] = -pt[1]; + grad[2] = 2 * quad->focal; +} + +static INLINE double inject_same_sign(double x, double src) { + const int64_t sign = (*(int64_t *) (&src)) & 0x8000000000000000; + const int64_t value = (*(int64_t *) (&x)) & 0x7FFFFFFFFFFFFFFF; + const int64_t result = sign | value; + return *(double *) (&result); +} + +int +solve_second + (const double a, + const double b, + const double c, + double* dist) +{ + ASSERT(dist); + if (a != 0) { + /* standard case: 2nd degree */ + const double delta = b * b - 4 * a * c; + if (delta > 0) { + const double sqrt_delta = sqrt(delta); + /* precise formula */ + const double t1 = (-b - inject_same_sign(sqrt_delta, b)) / (2 * b); + const double t2 = c / (a * t1); + if (t1 >= 0) { + *dist = t1 < t2 ? t1: t2; + return 1; + } + else if (t2 >= 0) { + *dist = t2; + } + return (t2 >= 0); + } + else if (delta == 0) { + const double t = -b / (2 * a); + if (t >= 0) *dist = t; + return (t >= 0); + } + else { + return 0; + } + } + else if (b != 0) { + /* degenerated case: 1st degree only */ + const double t = -c / b; + if (t >= 0) *dist = t; + return (t >= 0); + } + /* fully degenerated case: cannot determine dist */ + return 0; +} + +int +quadric_plane_intersect_local + (const double org[3], + const double dir[3], + double pt[3], + double grad[3], + double* dist) +{ + /* Define z = 0 */ + const double a = 0; + const double b = dir[2]; + const double c = org[2]; + int sol = solve_second(a, b, c, dist); + if (!sol) return 0; + d3_add(pt, org, d3_muld(pt, dir, *dist)); + quadric_plane_gradient_local(grad); + return 1; +} + +int +quadric_parabol_intersect_local + (const double org[3], + const double dir[3], + const struct ssol_quadric_parabol* quad, + double pt[3], + double grad[3], + double* dist) +{ + /* Define x^2 + y^2 - 4 focal z = 0 */ + const double a = dir[0] * dir[0] + dir[1] * dir[1]; + const double b = + 2 * org[0] * dir[0] + 2 * org[1] * dir[1] - 4 * quad->focal * dir[2]; + const double c = org[0] * org[0] + org[1] * org[1] - 4 * quad->focal * org[2]; + int sol = solve_second(a, b, c, dist); + if (!sol) return 0; + d3_add(pt, org, d3_muld(pt, dir, *dist)); + quadric_parabol_gradient_local(grad, pt, quad); + return 1; +} + +int +quadric_parabolic_cylinder_intersect_local + (const double org[3], + const double dir[3], + const struct ssol_quadric_parabolic_cylinder* quad, + double pt[3], + double grad[3], + double* dist) +{ + /* Define y^2 - 4 focal z = 0 */ + const double a = dir[1] * dir[1]; + const double b = 2 * org[1] * dir[1] - 4 * quad->focal * dir[2]; + const double c = org[1] * org[1] - 4 * quad->focal * org[2]; + int sol = solve_second(a, b, c, dist); + if (!sol) return 0; + d3_add(pt, org, d3_muld(pt, dir, *dist)); + quadric_parabolic_cylinder_gradient_local(grad, pt, quad); + return 1; +} + +/******************************************************************************* * Exported ssol_shape functions ******************************************************************************/ res_T diff --git a/src/ssol_shape_c.h b/src/ssol_shape_c.h @@ -35,4 +35,45 @@ struct ssol_shape { ref_T ref; }; +void +quadric_plane_gradient_local(double grad[3]); + +void +quadric_parabol_gradient_local + (double grad[3], + const double pt[3], + const struct ssol_quadric_parabol* quad); + +void +quadric_parabolic_cylinder_gradient_local + (double grad[3], + const double pt[3], + const struct ssol_quadric_parabolic_cylinder* quad); + +int +quadric_plane_intersect_local + (const double org[3], + const double dir[3], + double pt[3], + double grad[3], + double* dist); + +int +quadric_parabol_intersect_local + (const double org[3], + const double dir[3], + const struct ssol_quadric_parabol* quad, + double pt[3], + double grad[3], + double* dist); + +int +quadric_parabolic_cylinder_intersect_local + (const double org[3], + const double dir[3], + const struct ssol_quadric_parabolic_cylinder* quad, + double pt[3], + double grad[3], + double* dist); + #endif /* SSOL_SHAPE_C_H */ diff --git a/src/ssol_solver.c b/src/ssol_solver.c @@ -40,25 +40,13 @@ static const char* END_TEXT[] = END_TEXT__; /******************************************************************************* * Helper functions ******************************************************************************/ -static INLINE int -is_instance_punched(const struct ssol_instance* instance) -{ - ASSERT(instance); - return instance->object->shape->type == SHAPE_PUNCHED; -} - -static INLINE const struct ssol_quadric* -get_quadric(const struct ssol_instance* instance) -{ - ASSERT(instance && is_instance_punched(instance)); - return &instance->object->shape->quadric; -} static FINLINE void solstice_trace_ray(struct realisation* rs) { float org[3], dir[3]; struct segment* seg = current_segment(rs); + f3_set_d3(org, seg->org); f3_set_d3(dir, seg->dir); S3D(scene_view_trace_ray @@ -113,7 +101,7 @@ set_sun_distributions(struct solver_data* data) break; default: res = RES_OK; - FATAL("Unreachable code \n"); + FATAL("Unreachable code\n"); } exit: @@ -446,7 +434,7 @@ set_output_pos_and_dir(struct realisation* rs) { if (fst_segment) { const struct ssol_sun* sun = rs->data.scene->sun; ASSERT(-1 <= rs->start.cos_sun && rs->start.cos_sun <= 0); - seg->weight = sun_get_dni(sun) + seg->weight = sun->dni * brdf_composite_sample (rs->data.brdfs, rs->data.rng, rs->start.sundir, rs->data.fragment.Ns, seg->dir) * -rs->start.cos_sun; @@ -473,20 +461,6 @@ propagate(struct realisation* rs) ASSERT(get_material_from_hit(rs->data.scene, &seg->hit)->type != MATERIAL_VIRTUAL); - /* offset the impact point and recompute normal if needed */ - ASSERT(rs->data.instance); - switch (rs->data.instance->object->shape->type) { - case SHAPE_MESH: - /* no postprocess needed */ - break; - case SHAPE_PUNCHED: - /* project the impact point on the quadric */ - FATAL("TODO\n"); - /* compute normal to quadric */ - break; - default: FATAL("Unreachable code\n"); break; - } - /* fill fragment from hit and loop */ d3_muld(seg->hit_pos, seg->dir, (double)seg->hit.distance); d3_add(seg->hit_pos, seg->org, seg->hit_pos); diff --git a/src/ssol_solver_c.h b/src/ssol_solver_c.h @@ -62,6 +62,7 @@ struct segment { struct s3d_hit hit; double org[3], dir[4]; double hit_pos[3]; + double hit_normal[3]; }; struct starting_point { diff --git a/src/ssol_sun_c.h b/src/ssol_sun_c.h @@ -48,11 +48,4 @@ struct ssol_sun { ref_T ref; }; -static INLINE double -sun_get_dni(const struct ssol_sun* sun) -{ - ASSERT(sun); - return sun->dni; -} - #endif /* SSOL_SUN_C_H */