From 2bb0740b68fdb62d45adc00204b3990ca42ade77 Mon Sep 17 00:00:00 2001 From: Jaron Kent-Dobias Date: Mon, 22 Aug 2016 10:11:14 -0400 Subject: started repo again without all the data files gunking the works --- src/fortune/geometry.c | 184 +++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 184 insertions(+) create mode 100644 src/fortune/geometry.c (limited to 'src/fortune/geometry.c') diff --git a/src/fortune/geometry.c b/src/fortune/geometry.c new file mode 100644 index 0000000..131cacf --- /dev/null +++ b/src/fortune/geometry.c @@ -0,0 +1,184 @@ +# +#include "defs.h" +#include + +void geominit() { + struct Edge e; + double sn; + + freeinit(&efl, sizeof e); + nvertices = 0; + nedges = 0; + sn = nsites + 4; + sqrt_nsites = sqrt(sn); + deltay = ymax - ymin; + deltax = xmax - xmin; +} + +struct Edge *bisect(s1, s2) struct Site *s1, *s2; +{ + double dx, dy, adx, ady; + struct Edge *newedge; + + newedge = (struct Edge *)getfree(&efl); + + newedge->reg[0] = s1; + newedge->reg[1] = s2; + ref(s1); + ref(s2); + newedge->ep[0] = (struct Site *)NULL; + newedge->ep[1] = (struct Site *)NULL; + + dx = s2->coord.x - s1->coord.x; + dy = s2->coord.y - s1->coord.y; + adx = dx > 0 ? dx : -dx; + ady = dy > 0 ? dy : -dy; + newedge->c = s1->coord.x * dx + s1->coord.y * dy + (dx * dx + dy * dy) * 0.5; + if (adx > ady) { + newedge->a = 1.0; + newedge->b = dy / dx; + newedge->c /= dx; + } else { + newedge->b = 1.0; + newedge->a = dx / dy; + newedge->c /= dy; + }; + + newedge->edgenbr = nedges; + out_bisector(newedge); + nedges += 1; + return (newedge); +} + +struct Site *intersect(el1, el2) struct Halfedge *el1, *el2; +{ + struct Edge *e1, *e2, *e; + struct Halfedge *el; + double d, xint, yint; + int right_of_site; + struct Site *v; + + e1 = el1->ELedge; + e2 = el2->ELedge; + if (e1 == (struct Edge *)NULL || e2 == (struct Edge *)NULL) + return ((struct Site *)NULL); + if (e1->reg[1] == e2->reg[1]) + return ((struct Site *)NULL); + + d = e1->a * e2->b - e1->b * e2->a; + /* printf("intersect: d=%g\n", d); */ + if (-1.0e-20 < d && d < 1.0e-20) { + return ((struct Site *)NULL); + }; + + xint = (e1->c * e2->b - e2->c * e1->b) / d; + yint = (e2->c * e1->a - e1->c * e2->a) / d; + + if ((e1->reg[1]->coord.y < e2->reg[1]->coord.y) || + (e1->reg[1]->coord.y == e2->reg[1]->coord.y && + e1->reg[1]->coord.x < e2->reg[1]->coord.x)) { + el = el1; + e = e1; + } else { + el = el2; + e = e2; + }; + right_of_site = xint >= e->reg[1]->coord.x; + if ((right_of_site && el->ELpm == le) || (!right_of_site && el->ELpm == re)) + return ((struct Site *)NULL); + + v = (struct Site *)getfree(&sfl); + v->refcnt = 0; + v->coord.x = xint; + v->coord.y = yint; + return (v); +} + +/* returns 1 if p is to right of halfedge e */ +int right_of(el, p) struct Halfedge *el; +struct Point *p; +{ + struct Edge *e; + struct Site *topsite; + int right_of_site, above, fast; + double dxp, dyp, dxs, t1, t2, t3, yl; + + e = el->ELedge; + topsite = e->reg[1]; + right_of_site = p->x > topsite->coord.x; + if (right_of_site && el->ELpm == le) + return (1); + if (!right_of_site && el->ELpm == re) + return (0); + + if (e->a == 1.0) { + dyp = p->y - topsite->coord.y; + dxp = p->x - topsite->coord.x; + fast = 0; + if ((!right_of_site && e->b < 0.0) | (right_of_site && e->b >= 0.0)) { + above = dyp >= e->b * dxp; + fast = above; + } else { + above = p->x + p->y * e->b > e->c; + if (e->b < 0.0) + above = !above; + if (!above) + fast = 1; + }; + if (!fast) { + dxs = topsite->coord.x - (e->reg[0])->coord.x; + above = e->b * (dxp * dxp - dyp * dyp) < + dxs * dyp * (1.0 + 2.0 * dxp / dxs + e->b * e->b); + if (e->b < 0.0) + above = !above; + }; + } else /*e->b==1.0 */ + { + yl = e->c - e->a * p->x; + t1 = p->y - yl; + t2 = p->x - topsite->coord.x; + t3 = yl - topsite->coord.y; + above = t1 * t1 > t2 * t2 + t3 * t3; + }; + return (el->ELpm == le ? above : !above); +} + +void endpoint(e, lr, s) struct Edge *e; +int lr; +struct Site *s; +{ + e->ep[lr] = s; + ref(s); + if (e->ep[re - lr] == (struct Site *)NULL) { + } else { + out_ep(e); + deref(e->reg[le]); + deref(e->reg[re]); + makefree((struct Freenode *)e, &efl); + } +} + +double dist(s, t) struct Site *s, *t; +{ + double dx, dy; + dx = s->coord.x - t->coord.x; + dy = s->coord.y - t->coord.y; + return (sqrt(dx * dx + dy * dy)); +} + +void makevertex(v) struct Site *v; +{ + v->sitenbr = nvertices; + nvertices += 1; + out_vertex(v); +} + +void deref(v) struct Site *v; +{ + v->refcnt -= 1; + if (v->refcnt == 0) + makefree((struct Freenode *)v, &sfl); +} + +void ref(v) struct Site *v; +{ v->refcnt += 1; } -- cgit v1.2.3-70-g09d2