diff --git a/libs/util/mathlib.c b/libs/util/mathlib.c index 1fc03508b..f4bf46dba 100644 --- a/libs/util/mathlib.c +++ b/libs/util/mathlib.c @@ -35,6 +35,8 @@ # include #endif +#include "qfalloca.h" + #include #define IMPLEMENT_R_Cull @@ -42,6 +44,7 @@ #include "QF/mathlib.h" #include "QF/qtypes.h" +#include "QF/set.h" #include "QF/sys.h" VISIBLE int nanmask = 255 << 23; @@ -1348,12 +1351,13 @@ test_support_points(const vec_t **points, int *num_points, const vec3_t center) break; case 2: VectorSubtract (points[1], points[0], v); - VectorSubtract (center, points[0], d); + VectorAdd (points[0], points[1], d); + VectorScale (d, 0.5, d); + VectorSubtract (center, d, d); CrossProduct (v, d, n); nn = DotProduct (n, n); - dd = DotProduct (d, d); vv = DotProduct (v, v); - in_affine = nn < 1e-5 * vv * dd; + in_affine = nn < 1e-5 * vv * vv; break; case 3: VectorSubtract (points[1], points[0], a); @@ -1403,12 +1407,14 @@ test_support_points(const vec_t **points, int *num_points, const vec3_t center) sphere_t SmallestEnclosingBall (const vec3_t points[], int num_points) { + set_t was_support = SET_STATIC_INIT (num_points, alloca); sphere_t sphere; const vec_t *best; const vec_t *support[4]; int num_support; vec_t dist, best_dist; int i; + int best_i = 0; int iters = 0; if (num_points < 1) { @@ -1419,6 +1425,7 @@ SmallestEnclosingBall (const vec3_t points[], int num_points) for (i = 0; i < 4; i++) support[i] = 0; + set_empty (&was_support); VectorCopy (points[0], sphere.center); best_dist = dist = 0; @@ -1427,55 +1434,51 @@ SmallestEnclosingBall (const vec3_t points[], int num_points) dist = VectorDistance_fast (points[i], sphere.center); if (dist > best_dist) { best_dist = dist; + best_i = i; best = points[i]; } } num_support = 1; support[0] = best; sphere.radius = best_dist; // note: radius squared until the end + set_add (&was_support, best_i); while (!test_support_points (support, &num_support, sphere.center)) { - vec3_t affine; - vec3_t center_to_affine, center_to_point; - vec_t affine_dist, point_proj, point_dist, bound; - vec_t scale = 1; + vec3_t affine, v, p, r; + vec_t x, best_x = 0, rr, pv, pp; int i; - if (iters++ > 10) + if (iters++ > 2 * num_points) Sys_Error ("stuck SEB"); closest_affine_point (support, num_support, sphere.center, affine); - VectorSubtract (affine, sphere.center, center_to_affine); - affine_dist = DotProduct (center_to_affine, center_to_affine); - if (affine_dist < sphere.radius * 1e-5) { - // It's possible test_support_points failed due to precision - // issues - break; - } + VectorSubtract (support[0], affine, r); + rr = DotProduct (r, r); + VectorSubtract (sphere.center, affine, v); best = 0; for (i = 0; i < num_points; i++) { - if (points[i] == support[0] || points[i] == support[1] - || points[i] == support[2]) + if (SET_TEST_MEMBER (&was_support, i)) { continue; - VectorSubtract (points[i], sphere.center, center_to_point); - point_proj = DotProduct (center_to_affine, center_to_point); - if (affine_dist - point_proj <= 0 - || ((affine_dist - point_proj) * (affine_dist - point_proj) - < 1e-6 * sphere.radius * affine_dist)) + } + VectorSubtract (points[i], affine, p); + pp = DotProduct (p, p); + pv = DotProduct (p, v); + if (pp <= rr || pv <= 0 || pv * pv < 1e-6 * rr) { continue; - point_dist = DotProduct (center_to_point, center_to_point); - bound = sphere.radius - point_dist; - bound /= 2 * (affine_dist - point_proj); - if (bound < scale) { + } + x = (pp - rr) / (2 * pv); + if (x > best_x) { best = points[i]; - scale = bound; + best_i = i; + best_x = x; } } - VectorMultAdd (sphere.center, scale, center_to_affine, sphere.center); + VectorMultAdd (affine, best_x, v, sphere.center); sphere.radius = VectorDistance_fast (sphere.center, support[0]); if (best) { support[num_support++] = best; + set_add (&was_support, best_i); } } best_dist = 0; diff --git a/libs/util/simd.c b/libs/util/simd.c index 07670219b..5ad1ab3e0 100644 --- a/libs/util/simd.c +++ b/libs/util/simd.c @@ -28,6 +28,8 @@ # include "config.h" #endif +#include "qfalloca.h" + #include #define IMPLEMENT_VEC4F_Funcs @@ -37,6 +39,7 @@ #include "QF/simd/vec4d.h" #include "QF/simd/vec4f.h" #include "QF/simd/mat4f.h" +#include "QF/set.h" #include "QF/sys.h" vec4f_t @@ -192,12 +195,11 @@ test_support_points(const vec4f_t **points, int *num_points, vec4f_t center) break; case 2: v = *points[1] - *points[0]; - d = center - *points[0]; + d = center - (*points[0] + *points[1]) / 2; n = crossf (v, d); nn = dotf (n, n)[0]; - dd = dotf (d, d)[0]; vv = dotf (v, v)[0]; - in_affine = nn < 1e-5 * vv * dd; + in_affine = nn <= 1e-5 * vv * vv; break; case 3: a = *points[1] - *points[0]; @@ -247,6 +249,7 @@ test_support_points(const vec4f_t **points, int *num_points, vec4f_t center) vspheref_t SmallestEnclosingBall_vf (const vec4f_t *points, int num_points) { + set_t was_support = SET_STATIC_INIT (num_points, alloca); vspheref_t sphere = {}; vec4f_t center = {}; const vec4f_t *best; @@ -262,6 +265,7 @@ SmallestEnclosingBall_vf (const vec4f_t *points, int num_points) for (i = 0; i < 4; i++) { support[i] = 0; } + set_empty (&was_support); vec4f_t dist = {}; float best_dist = 0; @@ -278,50 +282,51 @@ SmallestEnclosingBall_vf (const vec4f_t *points, int num_points) num_support = 1; support[0] = best; sphere.radius = best_dist; // note: radius squared until the end + set_add (&was_support, best - points); while (!test_support_points (support, &num_support, center)) { - vec4f_t affine; - vec4f_t center_to_affine, center_to_point; - float affine_dist, point_proj, point_dist, bound; - float scale = 1; + vec4f_t affine, r, rr; + vec4f_t v, p, pv, pp, x; + vec4f_t best_x = { }; int i; - if (iters++ > 10) + //Sys_Printf ("%d: "VEC4F_FMT", %.9g, %d\n", iters, VEC4_EXP (center), sqrt(sphere.radius), num_support); + if (iters++ > 2 * num_points) { + //for (i = 0; i < num_points; i++) { + // Sys_Printf (VEC4F_FMT",\n", VEC4_EXP (points[i])); + //} Sys_Error ("stuck SEB"); + } affine = closest_affine_point (support, num_support, center); - center_to_affine = affine - center; - affine_dist = dotf (center_to_affine, center_to_affine)[0]; - if (affine_dist < sphere.radius * 1e-5) { - // It's possible test_support_points failed due to precision - // issues - break; - } + r = *support[0] - affine; //FIXME better choice + rr = dotf (r, r); + v = center - affine; best = 0; for (i = 0; i < num_points; i++) { - if (&points[i] == support[0] || &points[i] == support[1] - || &points[i] == support[2]) + if (SET_TEST_MEMBER (&was_support, i)) { continue; - center_to_point = points[i] - center; - point_proj = dotf (center_to_affine, center_to_point)[0]; - if (affine_dist - point_proj <= 0 - || ((affine_dist - point_proj) * (affine_dist - point_proj) - < 1e-6 * sphere.radius * affine_dist)) + } + p = points[i] - affine; + pp = dotf (p, p); + pv = dotf (p, v); + if (pp[0] <= rr[0] || pv[0] <= 0 + || (pv[0] * pv[0]) < 1e-6 * rr[0]) { continue; - point_dist = dotf (center_to_point, center_to_point)[0]; - bound = sphere.radius - point_dist; - bound /= 2 * (affine_dist - point_proj); - if (bound < scale) { + } + x = (pp - rr) / (2 * pv); + if (x[0] > best_x[0]) { best = &points[i]; - scale = bound; + best_x = x; } } - center = center + scale * center_to_affine; + center = affine + best_x * v; dist = center - *support[0]; sphere.radius = dotf (dist, dist)[0]; if (best) { support[num_support++] = best; + set_add (&was_support, best - points); } } best_dist = 0; diff --git a/libs/util/test/test-seb.c b/libs/util/test/test-seb.c index 01a4c6b25..a901e4f70 100644 --- a/libs/util/test/test-seb.c +++ b/libs/util/test/test-seb.c @@ -6,6 +6,8 @@ #include "QF/mersenne.h" #include "QF/sys.h" +#define SIZEOF(x) (sizeof(x) / sizeof(x[0])) + const vec3_t points[] = { {-1, -1, 1}, { 1, 1, 1}, @@ -24,9 +26,155 @@ const vec3_t points[] = { // due to qfbsp not culling the portal for being too small. const vec3_t tears_triangle[] = { {2201.82007, -1262, -713.450012}, - {2201.8501, -1262, -713.593994}, + {2201.8501, -1262, -713.593994}, {2201.84009, -1262, -713.445007}, }; +const vec3_t tears_triangle2[] = { + {-519.989014, 1111.94995, -2592.05005}, + {-520.002991, 1112, -2592.02002}, + {-520.005005, 1112, -2592.04004}, +}; +const vec3_t tears_triangle3[] = { + {284, 624, 2959.76001}, + {284, 624, 2960}, + {284.234985, 624.471008, 2959.53003}, +}; + +const vec3_t tears_quad[] = { + {1792, 1155.53003, 240}, + {1792, 1030.32996, 240}, + {1792.01001, 1030.32996, 240}, + {1792.01001, 1155.53003, 240}, +}; + +const vec3_t tears_quad2[] = { + {304, 3248, 3208}, + {-61.2307014, 3248, 3208}, + {-61.2307014, 3248, 3192}, + {304, 3248, 3192}, +}; + +const vec3_t tears_cluster[] = { + {551.638, -439.277008, 2804.63989}, + {551, -438, 2806.02002}, + {555.999023, -447.998993, 2805.1499}, + {555.999023, -447.998993, 2804.32007}, + {551.638, -439.277008, 2804.63989}, + {555.999023, -447.998993, 2804.32007}, + {555.999023, -447.998993, 2795.19995}, + {550.97699, -437.954987, 2806.02002}, + {551, -438, 2806.02002}, + {555.999023, -447.998993, 2795.19995}, + {555.999023, -447.998993, 2794.6499}, + {550.978027, -437.954987, 2806.02002}, + {545.945984, -442.945007, 2804.78003}, + {546, -442.998993, 2804.77002}, + {546, -443, 2804.69995}, + {546, -443, 2804.78003}, + {545.945007, -442.945007, 2804.78003}, + {556, -448, 2794.6499}, + {556, -448, 2805.1499}, + {546, -443, 2804.78003}, + {546, -443, 2804.69995}, + {550.97699, -437.954987, 2806.02002}, + {556, -448, 2794.6499}, + {546.000977, -443, 2804.69995}, + {545.945007, -442.945007, 2804.78003}, +}; + +const vec3_t tears_cluster2[] = { + {-2160, -192, -704.242004}, + {-2160, -192, -662}, + {-2160, -196, -658}, + {-2160, -220, -646}, + {-2160, -222.332001, -645.416992}, + {-2160, -320, -752}, + {-2160, -192, -752}, + {-2160, -192, -704.242004}, + {-2160, -222.332001, -645.416992}, + {-2160, -244, -640}, + {-2160, -268, -640}, + {-2160, -292, -646}, + {-2160, -316, -658}, + {-2160, -320, -662}, + {-2160, -196, -814}, + {-2160, -192, -810}, + {-2160, -192, -752}, + {-2160, -320, -752}, + {-2160, -320, -810}, + {-2160, -316, -814}, + {-2160, -312, -816}, + {-2160, -200, -816}, + {-2174, -334, -676}, + {-2186, -346, -700}, + {-2192, -352, -724}, + {-2192, -352, -748}, + {-2186, -346, -772}, + {-2174, -334, -796}, + {-2160, -319.998993, -810}, + {-2160, -319.998993, -661.999023}, + {-2174, -178, -796}, + {-2186, -166, -772}, + {-2192, -160, -748}, + {-2192, -160, -724}, + {-2186, -166, -700}, + {-2174, -178, -676}, + {-2160, -192, -662}, + {-2160, -192, -810}, + {-2224, -320, -810}, + {-2224, -330, -800}, + {-2224, -320, -800}, + {-2224, -200, -816}, + {-2224, -312, -816}, + {-2224, -316, -814}, + {-2224, -320, -810}, + {-2224, -320, -800}, + {-2224, -192, -800}, + {-2224, -192, -810}, + {-2224, -196, -814}, + {-2224, -192, -800}, + {-2224, -182, -800}, + {-2224, -192, -810}, + {-2224, -330, -800}, + {-2224, -334, -796}, + {-2224, -178, -796}, + {-2224, -182, -800}, + {-2224, -334, -796}, + {-2224, -346, -772}, + {-2224, -351, -752}, + {-2224, -161, -752}, + {-2224, -166, -772}, + {-2224, -178, -796}, + {-2224, -351, -752}, + {-2224, -352, -748}, + {-2224, -352, -724}, + {-2224, -351, -720}, + {-2224, -161, -720}, + {-2224, -160, -724}, + {-2224, -160, -748}, + {-2224, -161, -752}, + {-2224, -351, -720}, + {-2224, -346, -700}, + {-2224, -334, -676}, + {-2224, -316, -658}, + {-2224, -292, -646}, + {-2224, -268, -640}, + {-2224, -244, -640}, + {-2224, -220, -646}, + {-2224, -196, -658}, + {-2224, -178, -676}, + {-2224, -166, -700}, + {-2224, -161, -720}, +}; + +// random numbers dug this up. needed more than 6 iterations +const vec3_t cloud_points[] = { + {1.70695281, 2.60889673, 2.86233163}, + {2.11825109, 2.59883952, 1.48193693}, + {2.76728868, 1.57646465, 2.39769602}, + {1.88795376, 2.73980999, 1.58031297}, + {2.79736757, 1.89582968, 1.80514407}, +}; struct { const vec3_t *points; @@ -38,9 +186,24 @@ struct { {points, 2, {{ 0, 0, 1}, 1.41421356}}, {points, 3, {{-0.333333343, 0.333333343, 0.333333343}, 1.63299322}}, {points, 4, {{0, 0, 0}, 1.73205081}}, - {tears_triangle, 3, {{2201.84521, -1262, -713.519531}, 0.0747000724}}, + {tears_triangle, SIZEOF (tears_triangle), + {{2201.84521, -1262, -713.519531}, 0.0747000724}},//FIXME numbers? + {tears_cluster, SIZEOF (tears_cluster), + {{552.678101, -443.460876, 2800.40405}, 8.04672813}},//FIXME numbers? + {tears_quad, SIZEOF (tears_quad), + {{1792.005, 1092.92999, 240}, 62.6000302}}, + {tears_triangle2, SIZEOF (tears_triangle2), + {{-519.995972, 1111.97498, -2592.03516}, 0.0300767105}},//FIXME numbers? + {tears_quad2, SIZEOF (tears_quad2), + {{121.384659, 3248, 3200}, 182.790497}}, + {tears_triangle3, SIZEOF (tears_triangle3), + {{284.117493, 624.235535, 2959.76489}, 0.352925777}}, + {tears_cluster2, SIZEOF (tears_cluster2), + {{-2192, -256.000031, -736}, 103.479492}},//FIXME numbers? + {cloud_points, SIZEOF (cloud_points), + {{2.20056152, 2.23369908, 2.2332375}, 0.88327992}}, }; -#define num_tests (sizeof (tests) / sizeof (tests[0])) +#define num_tests SIZEOF (tests) static inline float rnd (mtstate_t *mt) diff --git a/libs/util/test/test-sebvf.c b/libs/util/test/test-sebvf.c index ef811e364..9b0d93a7f 100644 --- a/libs/util/test/test-sebvf.c +++ b/libs/util/test/test-sebvf.c @@ -7,6 +7,8 @@ #include "QF/sys.h" #include "QF/simd/vec4f.h" +#define SIZEOF(x) (sizeof(x) / sizeof(x[0])) + const vec4f_t points[] = { {-1, -1, 1, 1}, { 1, 1, 1, 1}, @@ -25,23 +27,186 @@ const vec4f_t points[] = { // due to qfbsp not culling the portal for being too small. const vec4f_t tears_triangle[] = { {2201.82007, -1262, -713.450012, 1}, - {2201.8501, -1262, -713.593994, 1}, + {2201.8501, -1262, -713.593994, 1}, {2201.84009, -1262, -713.445007, 1}, }; +const vec4f_t tears_triangle2[] = { + {-519.989014, 1111.94995, -2592.05005, 1}, + {-520.002991, 1112, -2592.02002, 1}, + {-520.005005, 1112, -2592.04004, 1}, +}; +const vec4f_t tears_triangle3[] = { + {284, 624, 2959.76001, 1}, + {284, 624, 2960, 1}, + {284.234985, 624.471008, 2959.53003, 1}, +}; + +const vec4f_t tears_quad[] = { + {1792, 1155.53003, 240, 1}, + {1792, 1030.32996, 240, 1}, + {1792.01001, 1030.32996, 240, 1}, + {1792.01001, 1155.53003, 240, 1}, +}; + +const vec4f_t tears_quad2[] = { + {304, 3248, 3208, 1}, + {-61.2307014, 3248, 3208, 1}, + {-61.2307014, 3248, 3192, 1}, + {304, 3248, 3192, 1}, +}; + +const vec4f_t tears_cluster[] = { + {551.638, -439.277008, 2804.63989, 1}, + {551, -438, 2806.02002, 1}, + {555.999023, -447.998993, 2805.1499, 1}, + {555.999023, -447.998993, 2804.32007, 1}, + {551.638, -439.277008, 2804.63989, 1}, + {555.999023, -447.998993, 2804.32007, 1}, + {555.999023, -447.998993, 2795.19995, 1}, + {550.97699, -437.954987, 2806.02002, 1}, + {551, -438, 2806.02002, 1}, + {555.999023, -447.998993, 2795.19995, 1}, + {555.999023, -447.998993, 2794.6499, 1}, + {550.978027, -437.954987, 2806.02002, 1}, + {545.945984, -442.945007, 2804.78003, 1}, + {546, -442.998993, 2804.77002, 1}, + {546, -443, 2804.69995, 1}, + {546, -443, 2804.78003, 1}, + {545.945007, -442.945007, 2804.78003, 1}, + {556, -448, 2794.6499, 1}, + {556, -448, 2805.1499, 1}, + {546, -443, 2804.78003, 1}, + {546, -443, 2804.69995, 1}, + {550.97699, -437.954987, 2806.02002, 1}, + {556, -448, 2794.6499, 1}, + {546.000977, -443, 2804.69995, 1}, + {545.945007, -442.945007, 2804.78003, 1}, +}; + +const vec4f_t tears_cluster2[] = { + {-2160, -192, -704.242004, 1}, + {-2160, -192, -662, 1}, + {-2160, -196, -658, 1}, + {-2160, -220, -646, 1}, + {-2160, -222.332001, -645.416992, 1}, + {-2160, -320, -752, 1}, + {-2160, -192, -752, 1}, + {-2160, -192, -704.242004, 1}, + {-2160, -222.332001, -645.416992, 1}, + {-2160, -244, -640, 1}, + {-2160, -268, -640, 1}, + {-2160, -292, -646, 1}, + {-2160, -316, -658, 1}, + {-2160, -320, -662, 1}, + {-2160, -196, -814, 1}, + {-2160, -192, -810, 1}, + {-2160, -192, -752, 1}, + {-2160, -320, -752, 1}, + {-2160, -320, -810, 1}, + {-2160, -316, -814, 1}, + {-2160, -312, -816, 1}, + {-2160, -200, -816, 1}, + {-2174, -334, -676, 1}, + {-2186, -346, -700, 1}, + {-2192, -352, -724, 1}, + {-2192, -352, -748, 1}, + {-2186, -346, -772, 1}, + {-2174, -334, -796, 1}, + {-2160, -319.998993, -810, 1}, + {-2160, -319.998993, -661.999023, 1}, + {-2174, -178, -796, 1}, + {-2186, -166, -772, 1}, + {-2192, -160, -748, 1}, + {-2192, -160, -724, 1}, + {-2186, -166, -700, 1}, + {-2174, -178, -676, 1}, + {-2160, -192, -662, 1}, + {-2160, -192, -810, 1}, + {-2224, -320, -810, 1}, + {-2224, -330, -800, 1}, + {-2224, -320, -800, 1}, + {-2224, -200, -816, 1}, + {-2224, -312, -816, 1}, + {-2224, -316, -814, 1}, + {-2224, -320, -810, 1}, + {-2224, -320, -800, 1}, + {-2224, -192, -800, 1}, + {-2224, -192, -810, 1}, + {-2224, -196, -814, 1}, + {-2224, -192, -800, 1}, + {-2224, -182, -800, 1}, + {-2224, -192, -810, 1}, + {-2224, -330, -800, 1}, + {-2224, -334, -796, 1}, + {-2224, -178, -796, 1}, + {-2224, -182, -800, 1}, + {-2224, -334, -796, 1}, + {-2224, -346, -772, 1}, + {-2224, -351, -752, 1}, + {-2224, -161, -752, 1}, + {-2224, -166, -772, 1}, + {-2224, -178, -796, 1}, + {-2224, -351, -752, 1}, + {-2224, -352, -748, 1}, + {-2224, -352, -724, 1}, + {-2224, -351, -720, 1}, + {-2224, -161, -720, 1}, + {-2224, -160, -724, 1}, + {-2224, -160, -748, 1}, + {-2224, -161, -752, 1}, + {-2224, -351, -720, 1}, + {-2224, -346, -700, 1}, + {-2224, -334, -676, 1}, + {-2224, -316, -658, 1}, + {-2224, -292, -646, 1}, + {-2224, -268, -640, 1}, + {-2224, -244, -640, 1}, + {-2224, -220, -646, 1}, + {-2224, -196, -658, 1}, + {-2224, -178, -676, 1}, + {-2224, -166, -700, 1}, + {-2224, -161, -720, 1}, +}; + +// random numbers dug this up. needed more than 6 iterations +const vec4f_t cloud_points[] = { + {1.70695281, 2.60889673, 2.86233163, 1}, + {2.11825109, 2.59883952, 1.48193693, 1}, + {2.76728868, 1.57646465, 2.39769602, 1}, + {1.88795376, 2.73980999, 1.58031297, 1}, + {2.79736757, 1.89582968, 1.80514407, 1}, +}; struct { const vec4f_t *points; int num_points; vspheref_t expect; } tests[] = { +#if 0 {0, 0, {{ 0, 0, 0, 1}, 0}}, {points, 1, {{-1, -1, 1, 1}, 0}}, {points, 2, {{ 0, 0, 1, 1}, 1.41421356}}, {points, 3, {{-0.333333343, 0.333333343, 0.333333343, 1}, 1.63299322}}, {points, 4, {{0, 0, 0, 1}, 1.73205081}}, - {tears_triangle, 3, {{2201.84521, -1262, -713.519531}, 0.0747000724}}, + {tears_triangle, SIZEOF (tears_triangle), + {{2201.84521, -1262, -713.519531}, 0.0747000724}},//FIXME numbers? + {tears_cluster, SIZEOF (tears_cluster), + {{552.678101, -443.460876, 2800.40405}, 8.04672813}},//FIXME numbers? + {tears_quad, SIZEOF (tears_quad), + {{1792.005, 1092.92999, 240}, 62.6000302}}, + {tears_triangle2, SIZEOF (tears_triangle2), + {{-519.995972, 1111.97498, -2592.03516}, 0.0300767105}},//FIXME numbers? + {tears_quad2, SIZEOF (tears_quad2), + {{121.384659, 3248, 3200}, 182.790497}}, + {tears_triangle3, SIZEOF (tears_triangle3), + {{284.117493, 624.235535, 2959.76489}, 0.352925777}}, + {tears_cluster2, SIZEOF (tears_cluster2), + {{-2192, -256.000031, -736}, 103.479492}},//FIXME numbers? +#endif + {cloud_points, SIZEOF (cloud_points), + {{2.20056152, 2.23369908, 2.2332375}, 0.88327992}}, }; -#define num_tests (sizeof (tests) / sizeof (tests[0])) +#define num_tests SIZEOF (tests) static inline float rnd (mtstate_t *mt) @@ -92,6 +257,7 @@ main (int argc, const char **argv) for (j = 0; j < 5; j++) { VectorSet (rnd (&mt), rnd (&mt), rnd (&mt), cloud[j]); + cloud[j][3] = 1; } seb = SmallestEnclosingBall_vf (cloud, 5); r2 = seb.radius * seb.radius;