solstice-solver

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

commit 8fa91476e1c324465cfdcf0f5730d14fa598d7a7
parent 30aed1dc68b07f7d49db222954e40c2a6793a510
Author: Christophe Coustet <christophe.coustet@meso-star.com>
Date:   Wed,  7 Sep 2016 14:25:55 +0200

Fix ray/quadric intersection.

Previous code returned the neerest intersection point from origin.
Now we select the closest solution to the raytracing result.

Some changes in arg order in various functions, though.

Diffstat:
Msrc/ssol_instance_c.h | 8+-------
Msrc/ssol_scene.c | 18+++++++++++++-----
Msrc/ssol_shape.c | 67++++++++++++++++++++++++++++++++++++++++++-------------------------
Msrc/ssol_shape_c.h | 23++++++++++-------------
Msrc/ssol_solver_c.h | 1+
5 files changed, 67 insertions(+), 50 deletions(-)

diff --git a/src/ssol_instance_c.h b/src/ssol_instance_c.h @@ -19,6 +19,7 @@ #include <rsys/list.h> #include <rsys/ref_count.h> #include <rsys/str.h> + struct ssol_instance { struct ssol_object* object; /* Instantiated object */ struct s3d_shape* shape_rt; /* Instantiated Star-3D shape to ray-trace */ @@ -32,11 +33,4 @@ struct ssol_instance { ref_T ref; }; -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 @@ -294,13 +294,12 @@ hit_filter_function inst = *htable_instance_find(&rs->data.scene->instances_rt, &hit->prim.inst_id); - if (is_instance_punched(inst)) { + if (inst->object->shape->type == SHAPE_PUNCHED) { /* 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); + seg->org, seg->dir, seg->hit_pos, seg->hit_normal, &seg->dist); if (!ok) return 1; /* invalid impact */ break; } @@ -308,7 +307,7 @@ hit_filter_function 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); + quad, seg->org, seg->dir, hit->distance, seg->hit_pos, seg->hit_normal, &seg->dist); if (!ok) return 1; /* invalid impact */ break; } @@ -316,13 +315,22 @@ hit_filter_function 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); + quad, seg->org, seg->dir, hit->distance, seg->hit_pos, seg->hit_normal, &seg->dist); if (!ok) return 1; /* invalid impact */ break; } default: FATAL("Unreachable code\n"); break; } } + else { + ASSERT(inst->object->shape->type == SHAPE_MESH); + /* just copy raytracing results to segment */ + float pos[3]; + f3_add(pos, org, f3_mulf(pos, dir, hit->distance)); + d3_set_f3(seg->hit_pos, pos); + d3_set_f3(seg->hit_normal, hit->normal); + seg->dist = hit->distance; + } front_face = d3_dot(seg->hit_normal, seg->dir) < 0; diff --git a/src/ssol_shape.c b/src/ssol_shape.c @@ -537,8 +537,16 @@ inject_same_sign(const double x, const double src) return ucast.d; } +/* solve a 2nd degree equation + hint is used to select among the 2 solutions (if applies) + the selected solution is then the closest to hint positive value */ static int -quadric_solve_second(const double a, const double b, const double c, double* dist) +quadric_solve_second + (const double a, + const double b, + const double c, + const double hint, + double* dist) { ASSERT(dist); if (a != 0) { @@ -549,25 +557,32 @@ quadric_solve_second(const double a, const double b, const double c, double* dis /* 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; + if (t1 < 0 && t2 < 0) return 0; /* no positive solution */ + if (t1 < 0) { + *dist = t2; /* t2 is the only positive solution */ + return 1; + } + if (t2 < 0) { + *dist = t1; /* t1 is the only positive solution */ return 1; - } else if (t2 >= 0) { - *dist = t2; } - return (t2 >= 0); + /* both t1 and t2 are positive: choose the closest value to hint */ + *dist = fabs(t1 - hint) < fabs(t2 - hint) ? t1 : t2; + return 1; } else if (delta == 0) { const double t = -b / (2 * a); - if (t >= 0) *dist = t; - return (t >= 0); + if (t < 0) return 0; /* no positive solution */ + *dist = t; + return 1; } 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); + if (t < 0) return 0; /* no positive solution */ + *dist = t; + return 1; } /* fully degenerated case: cannot determine dist */ return 0; @@ -600,9 +615,9 @@ quadric_plane_gradient_local(double grad[3]) void quadric_parabol_gradient_local - (double grad[3], + (const struct ssol_quadric_parabol* quad, const double pt[3], - const struct ssol_quadric_parabol* quad) + double grad[3]) { grad[0] = -pt[0]; grad[1] = -pt[1]; @@ -611,9 +626,9 @@ quadric_parabol_gradient_local void quadric_parabolic_cylinder_gradient_local - (double grad[3], + (const struct ssol_quadric_parabolic_cylinder* quad, const double pt[3], - const struct ssol_quadric_parabolic_cylinder* quad) + double grad[3]) { grad[0] = 0; grad[1] = -pt[1]; @@ -628,11 +643,11 @@ quadric_plane_intersect_local double grad[3], double* dist) { - /* Define z = 0 */ + /* Define 0 z^2 + z + 0 = 0 */ const double a = 0; const double b = dir[2]; const double c = org[2]; - int sol = quadric_solve_second(a, b, c, dist); + int sol = quadric_solve_second(a, b, c, 0, dist); if (!sol) return 0; d3_add(pt, org, d3_muld(pt, dir, *dist)); quadric_plane_gradient_local(grad); @@ -641,30 +656,32 @@ quadric_plane_intersect_local int quadric_parabol_intersect_local - (const double org[3], + (const struct ssol_quadric_parabol* quad, + const double org[3], const double dir[3], - const struct ssol_quadric_parabol* quad, + const double hint, double pt[3], double grad[3], - double* dist) + double* dist) /* in/out: */ { /* 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]; - const int sol = quadric_solve_second(a, b, c, dist); + const int sol = quadric_solve_second(a, b, c, hint, dist); if (!sol) return 0; d3_add(pt, org, d3_muld(pt, dir, *dist)); - quadric_parabol_gradient_local(grad, pt, quad); + quadric_parabol_gradient_local(quad, pt, grad); return 1; } int quadric_parabolic_cylinder_intersect_local - (const double org[3], + (const struct ssol_quadric_parabolic_cylinder* quad, + const double org[3], const double dir[3], - const struct ssol_quadric_parabolic_cylinder* quad, + const double hint, double pt[3], double grad[3], double* dist) @@ -673,10 +690,10 @@ quadric_parabolic_cylinder_intersect_local 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]; - const int sol = quadric_solve_second(a, b, c, dist); + const int sol = quadric_solve_second(a, b, c, hint, dist); if (!sol) return 0; d3_add(pt, org, d3_muld(pt, dir, *dist)); - quadric_parabolic_cylinder_gradient_local(grad, pt, quad); + quadric_parabolic_cylinder_gradient_local(quad, pt, grad); return 1; } diff --git a/src/ssol_shape_c.h b/src/ssol_shape_c.h @@ -35,26 +35,21 @@ struct ssol_shape { ref_T ref; }; -/* - * FIXME: functions whose prefix is "quadric" should have the quadric data as - * first argument. - */ - extern LOCAL_SYM void quadric_plane_gradient_local (double grad[3]); extern LOCAL_SYM void quadric_parabol_gradient_local - (double grad[3], + (const struct ssol_quadric_parabol* quad, const double pt[3], - const struct ssol_quadric_parabol* quad); + double grad[3]); extern LOCAL_SYM void quadric_parabolic_cylinder_gradient_local - (double grad[3], + (const struct ssol_quadric_parabolic_cylinder* quad, const double pt[3], - const struct ssol_quadric_parabolic_cylinder* quad); + double grad[3]); extern LOCAL_SYM int quadric_plane_intersect_local @@ -66,18 +61,20 @@ quadric_plane_intersect_local extern LOCAL_SYM int quadric_parabol_intersect_local - (const double org[3], + (const struct ssol_quadric_parabol* quad, + const double org[3], const double dir[3], - const struct ssol_quadric_parabol* quad, + const double hint, double pt[3], double grad[3], double* dist); extern LOCAL_SYM int quadric_parabolic_cylinder_intersect_local - (const double org[3], + (const struct ssol_quadric_parabolic_cylinder* quad, + const double org[3], const double dir[3], - const struct ssol_quadric_parabolic_cylinder* quad, + const double hint, double pt[3], double grad[3], double* dist); diff --git a/src/ssol_solver_c.h b/src/ssol_solver_c.h @@ -63,6 +63,7 @@ struct segment { double org[3], dir[4]; double hit_pos[3]; double hit_normal[3]; + double dist; }; struct starting_point {