summaryrefslogtreecommitdiff
path: root/meme/slice.c
diff options
context:
space:
mode:
authorWerner Almesberger <werner@almesberger.net>2015-03-14 19:35:58 (GMT)
committerWerner Almesberger <werner@almesberger.net>2015-03-14 19:35:58 (GMT)
commit6b7491bd5298dd31c09c63bf08adbc1f7013345d (patch)
tree436ff9e7866b8da6f95db6108e2d9284f8c99a3a /meme/slice.c
parent21d686b4c66f725ddf76e9e5fa6a416e7028fe77 (diff)
downloadmisc-6b7491bd5298dd31c09c63bf08adbc1f7013345d.zip
misc-6b7491bd5298dd31c09c63bf08adbc1f7013345d.tar.gz
misc-6b7491bd5298dd31c09c63bf08adbc1f7013345d.tar.bz2
meme/: new option -n to select outline calculation based on needle shape
Diffstat (limited to 'meme/slice.c')
-rw-r--r--meme/slice.c80
1 files changed, 71 insertions, 9 deletions
diff --git a/meme/slice.c b/meme/slice.c
index 9bf2ce5..2ce8657 100644
--- a/meme/slice.c
+++ b/meme/slice.c
@@ -21,7 +21,8 @@
#include "slice.h"
-#define EPSILON 1e-6
+#define EPSILON 1e-6
+#define MAX_WARN 100
struct point {
@@ -38,7 +39,6 @@ static struct point *points = NULL;
static double z0, cx, cy; /* z = cx * x + cy * y + z0 */
-
static void slice_plane(void)
{
/* @@@ */
@@ -87,9 +87,9 @@ static bool solve(const struct vertex *va, const struct vertex *vb, double *res)
return 1;
}
+
static void calc(void (*fn)(double x, double y, double z, bool last),
- const struct vertex *va, const struct vertex *vb, double t,
- bool needle, bool last)
+ const struct vertex *va, const struct vertex *vb, double t, bool last)
{
double dx, dy, dz;
@@ -97,22 +97,84 @@ static void calc(void (*fn)(double x, double y, double z, bool last),
dy = vb->y - va->y;
dz = vb->z - va->z;
- /* @@@ needle correction */
-
fn((va->x + t * dx) / 1000, (va->y + t * dy) / 1000,
(va->z + t * dz) / 1000, last);
}
+/*
+ * Known bug: we don't take into account horizontal needle deflection.
+ * E.g., when slightly beyond a vertical drop, the needle may bend as
+ * it plunges, thus reporting a lower height than a perfectly rigid
+ * needle would.
+ *
+ * If this happens, e may become larger than d and we simply go back to
+ * using the last point before the edge.
+ */
+
+
+#define NEEDLE_TIP_R 40 /* 40 um */
+#define NEEDLE_SIDE 14 /* dy = 14 * dx */
+
+
+static void step(void (*fn)(double x, double y, double z, bool last),
+ const struct vertex *va, const struct vertex *vb, bool last)
+{
+ static unsigned warnings = 0;
+ double dx, dy, dz;
+ double d, e, t;
+
+ if (va->z > vb->z) {
+ step(fn, vb, va, last);
+ return;
+ }
+
+ dx = vb->x - va->x;
+ dy = vb->y - va->y;
+ dz = vb->z - va->z;
+
+ d = hypot(dx, dy);
+
+ if (dz < NEEDLE_TIP_R) {
+ /*
+ * From:
+ * e^2 + (r - z)^2 = r^2
+ */
+ e = sqrt(2 * NEEDLE_TIP_R * dz - dz * dz);
+ } else {
+ e = (dz - NEEDLE_TIP_R) / NEEDLE_SIDE + NEEDLE_TIP_R;
+ }
+
+ if (e > d) {
+ if (warnings < MAX_WARN) {
+ fprintf(stderr, "e > d: %f, %f; dz = %f\n", e, d, dz);
+ } else if (warnings == MAX_WARN) {
+ fprintf(stderr, "...\n");
+ warnings++;
+ }
+ fn(vb->x / 1000.0, vb->y / 1000.0, vb->z / 1000.0, last);
+ } else {
+ t = fabs(d) < EPSILON ? 0 : e / d;
+ fn((va->x + t * dx) / 1000, (va->y + t * dy) / 1000,
+ vb->z / 1000, last);
+ }
+}
+
+
static void segment(void (*fn)(double x, double y, double z, bool last),
const struct vertex *va, const struct vertex *vb, const struct vertex *vc,
bool needle)
{
double tb, tc;
- if (solve(va, vb, &tb) && solve(va, vc, &tc)) {
- calc(fn, va, vb, tb, needle, 0);
- calc(fn, va, vc, tc, needle, 1);
+ if (needle) {
+ step(fn, va, vb, 0);
+ step(fn, va, vc, 1);
+ } else {
+ if (solve(va, vb, &tb) && solve(va, vc, &tc)) {
+ calc(fn, va, vb, tb, 0);
+ calc(fn, va, vc, tc, 1);
+ }
}
}