g->edges = NULL;
g->dots = NULL;
g->num_faces = g->num_edges = g->num_dots = 0;
- g->middle_face = NULL;
g->refcount = 1;
g->lowest_x = g->lowest_y = g->highest_x = g->highest_y = 0;
return g;
* Returns the nearest edge, or NULL if no edge is reasonably
* near the position.
*
- * This algorithm is nice and generic, and doesn't depend on any particular
- * geometric layout of the grid:
- * Start at any dot (pick one next to middle_face).
- * Walk along a path by choosing, from all nearby dots, the one that is
- * nearest the target (x,y). Hopefully end up at the dot which is closest
- * to (x,y). Should work, as long as faces aren't too badly shaped.
- * Then examine each edge around this dot, and pick whichever one is
- * closest (perpendicular distance) to (x,y).
- * Using perpendicular distance is not quite right - the edge might be
- * "off to one side". So we insist that the triangle with (x,y) has
- * acute angles at the edge's dots.
+ * Just judging edges by perpendicular distance is not quite right -
+ * the edge might be "off to one side". So we insist that the triangle
+ * with (x,y) has acute angles at the edge's dots.
*
* edge1
* *---------*------
*/
grid_edge *grid_nearest_edge(grid *g, int x, int y)
{
- grid_dot *cur;
grid_edge *best_edge;
double best_distance = 0;
int i;
- cur = g->middle_face->dots[0];
-
- for (;;) {
- /* Target to beat */
- long dist = SQ((long)cur->x - (long)x) + SQ((long)cur->y - (long)y);
- /* Look for nearer dot - if found, store in 'new'. */
- grid_dot *new = cur;
- int i;
- /* Search all dots in all faces touching this dot. Some shapes
- * (such as in Cairo) don't quite work properly if we only search
- * the dot's immediate neighbours. */
- for (i = 0; i < cur->order; i++) {
- grid_face *f = cur->faces[i];
- int j;
- if (!f) continue;
- for (j = 0; j < f->order; j++) {
- long new_dist;
- grid_dot *d = f->dots[j];
- if (d == cur) continue;
- new_dist = SQ((long)d->x - (long)x) + SQ((long)d->y - (long)y);
- if (new_dist < dist) {
- new = d;
- break; /* found closer dot */
- }
- }
- if (new != cur)
- break; /* found closer dot */
- }
-
- if (new == cur) {
- /* Didn't find a closer dot among the neighbours of 'cur' */
- break;
- } else {
- cur = new;
- }
- }
- /* 'cur' is nearest dot, so find which of the dot's edges is closest. */
best_edge = NULL;
- for (i = 0; i < cur->order; i++) {
- grid_edge *e = cur->edges[i];
+ for (i = 0; i < g->num_edges; i++) {
+ grid_edge *e = &g->edges[i];
long e2; /* squared length of edge */
long a2, b2; /* squared lengths of other sides */
double dist;
}
printf("]\n");
}
- printf("Middle face: %d\n", (int)(g->middle_face - g->faces));
}
/* Show the derived grid information, computed by grid_make_consistent */
static void grid_print_derived(grid *g)
for (i = 0; i < face_size; i++)
new_face->dots[i] = NULL;
new_face->edges = NULL;
+ new_face->has_incentre = FALSE;
g->num_faces++;
}
/* Assumes dot list has enough space */
* Assumes g->dots has enough capacity allocated */
static grid_dot *grid_get_dot(grid *g, tree234 *dot_list, int x, int y)
{
- grid_dot test = {0, NULL, NULL, x, y};
- grid_dot *ret = find234(dot_list, &test, NULL);
+ grid_dot test, *ret;
+
+ test.order = 0;
+ test.edges = NULL;
+ test.faces = NULL;
+ test.x = x;
+ test.y = y;
+ ret = find234(dot_list, &test, NULL);
if (ret)
return ret;
last_face->dots[position] = d;
}
+/*
+ * Helper routines for grid_find_incentre.
+ */
+static int solve_2x2_matrix(double mx[4], double vin[2], double vout[2])
+{
+ double inv[4];
+ double det;
+ det = (mx[0]*mx[3] - mx[1]*mx[2]);
+ if (det == 0)
+ return FALSE;
+
+ inv[0] = mx[3] / det;
+ inv[1] = -mx[1] / det;
+ inv[2] = -mx[2] / det;
+ inv[3] = mx[0] / det;
+
+ vout[0] = inv[0]*vin[0] + inv[1]*vin[1];
+ vout[1] = inv[2]*vin[0] + inv[3]*vin[1];
+
+ return TRUE;
+}
+static int solve_3x3_matrix(double mx[9], double vin[3], double vout[3])
+{
+ double inv[9];
+ double det;
+
+ det = (mx[0]*mx[4]*mx[8] + mx[1]*mx[5]*mx[6] + mx[2]*mx[3]*mx[7] -
+ mx[0]*mx[5]*mx[7] - mx[1]*mx[3]*mx[8] - mx[2]*mx[4]*mx[6]);
+ if (det == 0)
+ return FALSE;
+
+ inv[0] = (mx[4]*mx[8] - mx[5]*mx[7]) / det;
+ inv[1] = (mx[2]*mx[7] - mx[1]*mx[8]) / det;
+ inv[2] = (mx[1]*mx[5] - mx[2]*mx[4]) / det;
+ inv[3] = (mx[5]*mx[6] - mx[3]*mx[8]) / det;
+ inv[4] = (mx[0]*mx[8] - mx[2]*mx[6]) / det;
+ inv[5] = (mx[2]*mx[3] - mx[0]*mx[5]) / det;
+ inv[6] = (mx[3]*mx[7] - mx[4]*mx[6]) / det;
+ inv[7] = (mx[1]*mx[6] - mx[0]*mx[7]) / det;
+ inv[8] = (mx[0]*mx[4] - mx[1]*mx[3]) / det;
+
+ vout[0] = inv[0]*vin[0] + inv[1]*vin[1] + inv[2]*vin[2];
+ vout[1] = inv[3]*vin[0] + inv[4]*vin[1] + inv[5]*vin[2];
+ vout[2] = inv[6]*vin[0] + inv[7]*vin[1] + inv[8]*vin[2];
+
+ return TRUE;
+}
+
+void grid_find_incentre(grid_face *f)
+{
+ double xbest, ybest, bestdist;
+ int i, j, k, m;
+ grid_dot *edgedot1[3], *edgedot2[3];
+ grid_dot *dots[3];
+ int nedges, ndots;
+
+ if (f->has_incentre)
+ return;
+
+ /*
+ * Find the point in the polygon with the maximum distance to any
+ * edge or corner.
+ *
+ * Such a point must exist which is in contact with at least three
+ * edges and/or vertices. (Proof: if it's only in contact with two
+ * edges and/or vertices, it can't even be at a _local_ maximum -
+ * any such circle can always be expanded in some direction.) So
+ * we iterate through all 3-subsets of the combined set of edges
+ * and vertices; for each subset we generate one or two candidate
+ * points that might be the incentre, and then we vet each one to
+ * see if it's inside the polygon and what its maximum radius is.
+ *
+ * (There's one case which this algorithm will get noticeably
+ * wrong, and that's when a continuum of equally good answers
+ * exists due to parallel edges. Consider a long thin rectangle,
+ * for instance, or a parallelogram. This algorithm will pick a
+ * point near one end, and choose the end arbitrarily; obviously a
+ * nicer point to choose would be in the centre. To fix this I
+ * would have to introduce a special-case system which detected
+ * parallel edges in advance, set aside all candidate points
+ * generated using both edges in a parallel pair, and generated
+ * some additional candidate points half way between them. Also,
+ * of course, I'd have to cope with rounding error making such a
+ * point look worse than one of its endpoints. So I haven't done
+ * this for the moment, and will cross it if necessary when I come
+ * to it.)
+ *
+ * We don't actually iterate literally over _edges_, in the sense
+ * of grid_edge structures. Instead, we fill in edgedot1[] and
+ * edgedot2[] with a pair of dots adjacent in the face's list of
+ * vertices. This ensures that we get the edges in consistent
+ * orientation, which we could not do from the grid structure
+ * alone. (A moment's consideration of an order-3 vertex should
+ * make it clear that if a notional arrow was written on each
+ * edge, _at least one_ of the three faces bordering that vertex
+ * would have to have the two arrows tip-to-tip or tail-to-tail
+ * rather than tip-to-tail.)
+ */
+ nedges = ndots = 0;
+ bestdist = 0;
+ xbest = ybest = 0;
+
+ for (i = 0; i+2 < 2*f->order; i++) {
+ if (i < f->order) {
+ edgedot1[nedges] = f->dots[i];
+ edgedot2[nedges++] = f->dots[(i+1)%f->order];
+ } else
+ dots[ndots++] = f->dots[i - f->order];
+
+ for (j = i+1; j+1 < 2*f->order; j++) {
+ if (j < f->order) {
+ edgedot1[nedges] = f->dots[j];
+ edgedot2[nedges++] = f->dots[(j+1)%f->order];
+ } else
+ dots[ndots++] = f->dots[j - f->order];
+
+ for (k = j+1; k < 2*f->order; k++) {
+ double cx[2], cy[2]; /* candidate positions */
+ int cn = 0; /* number of candidates */
+
+ if (k < f->order) {
+ edgedot1[nedges] = f->dots[k];
+ edgedot2[nedges++] = f->dots[(k+1)%f->order];
+ } else
+ dots[ndots++] = f->dots[k - f->order];
+
+ /*
+ * Find a point, or pair of points, equidistant from
+ * all the specified edges and/or vertices.
+ */
+ if (nedges == 3) {
+ /*
+ * Three edges. This is a linear matrix equation:
+ * each row of the matrix represents the fact that
+ * the point (x,y) we seek is at distance r from
+ * that edge, and we solve three of those
+ * simultaneously to obtain x,y,r. (We ignore r.)
+ */
+ double matrix[9], vector[3], vector2[3];
+ int m;
+
+ for (m = 0; m < 3; m++) {
+ int x1 = edgedot1[m]->x, x2 = edgedot2[m]->x;
+ int y1 = edgedot1[m]->y, y2 = edgedot2[m]->y;
+ int dx = x2-x1, dy = y2-y1;
+
+ /*
+ * ((x,y) - (x1,y1)) . (dy,-dx) = r |(dx,dy)|
+ *
+ * => x dy - y dx - r |(dx,dy)| = (x1 dy - y1 dx)
+ */
+ matrix[3*m+0] = dy;
+ matrix[3*m+1] = -dx;
+ matrix[3*m+2] = -sqrt((double)dx*dx+(double)dy*dy);
+ vector[m] = (double)x1*dy - (double)y1*dx;
+ }
+
+ if (solve_3x3_matrix(matrix, vector, vector2)) {
+ cx[cn] = vector2[0];
+ cy[cn] = vector2[1];
+ cn++;
+ }
+ } else if (nedges == 2) {
+ /*
+ * Two edges and a dot. This will end up in a
+ * quadratic equation.
+ *
+ * First, look at the two edges. Having our point
+ * be some distance r from both of them gives rise
+ * to a pair of linear equations in x,y,r of the
+ * form
+ *
+ * (x-x1) dy - (y-y1) dx = r sqrt(dx^2+dy^2)
+ *
+ * We eliminate r between those equations to give
+ * us a single linear equation in x,y describing
+ * the locus of points equidistant from both lines
+ * - i.e. the angle bisector.
+ *
+ * We then choose one of x,y to be a parameter t,
+ * and derive linear formulae for x,y,r in terms
+ * of t. This enables us to write down the
+ * circular equation (x-xd)^2+(y-yd)^2=r^2 as a
+ * quadratic in t; solving that and substituting
+ * in for x,y gives us two candidate points.
+ */
+ double eqs[2][4]; /* a,b,c,d : ax+by+cr=d */
+ double eq[3]; /* a,b,c: ax+by=c */
+ double xt[2], yt[2], rt[2]; /* a,b: {x,y,r}=at+b */
+ double q[3]; /* a,b,c: at^2+bt+c=0 */
+ double disc;
+
+ /* Find equations of the two input lines. */
+ for (m = 0; m < 2; m++) {
+ int x1 = edgedot1[m]->x, x2 = edgedot2[m]->x;
+ int y1 = edgedot1[m]->y, y2 = edgedot2[m]->y;
+ int dx = x2-x1, dy = y2-y1;
+
+ eqs[m][0] = dy;
+ eqs[m][1] = -dx;
+ eqs[m][2] = -sqrt(dx*dx+dy*dy);
+ eqs[m][3] = x1*dy - y1*dx;
+ }
+
+ /* Derive the angle bisector by eliminating r. */
+ eq[0] = eqs[0][0]*eqs[1][2] - eqs[1][0]*eqs[0][2];
+ eq[1] = eqs[0][1]*eqs[1][2] - eqs[1][1]*eqs[0][2];
+ eq[2] = eqs[0][3]*eqs[1][2] - eqs[1][3]*eqs[0][2];
+
+ /* Parametrise x and y in terms of some t. */
+ if (abs(eq[0]) < abs(eq[1])) {
+ /* Parameter is x. */
+ xt[0] = 1; xt[1] = 0;
+ yt[0] = -eq[0]/eq[1]; yt[1] = eq[2]/eq[1];
+ } else {
+ /* Parameter is y. */
+ yt[0] = 1; yt[1] = 0;
+ xt[0] = -eq[1]/eq[0]; xt[1] = eq[2]/eq[0];
+ }
+
+ /* Find a linear representation of r using eqs[0]. */
+ rt[0] = -(eqs[0][0]*xt[0] + eqs[0][1]*yt[0])/eqs[0][2];
+ rt[1] = (eqs[0][3] - eqs[0][0]*xt[1] -
+ eqs[0][1]*yt[1])/eqs[0][2];
+
+ /* Construct the quadratic equation. */
+ q[0] = -rt[0]*rt[0];
+ q[1] = -2*rt[0]*rt[1];
+ q[2] = -rt[1]*rt[1];
+ q[0] += xt[0]*xt[0];
+ q[1] += 2*xt[0]*(xt[1]-dots[0]->x);
+ q[2] += (xt[1]-dots[0]->x)*(xt[1]-dots[0]->x);
+ q[0] += yt[0]*yt[0];
+ q[1] += 2*yt[0]*(yt[1]-dots[0]->y);
+ q[2] += (yt[1]-dots[0]->y)*(yt[1]-dots[0]->y);
+
+ /* And solve it. */
+ disc = q[1]*q[1] - 4*q[0]*q[2];
+ if (disc >= 0) {
+ double t;
+
+ disc = sqrt(disc);
+
+ t = (-q[1] + disc) / (2*q[0]);
+ cx[cn] = xt[0]*t + xt[1];
+ cy[cn] = yt[0]*t + yt[1];
+ cn++;
+
+ t = (-q[1] - disc) / (2*q[0]);
+ cx[cn] = xt[0]*t + xt[1];
+ cy[cn] = yt[0]*t + yt[1];
+ cn++;
+ }
+ } else if (nedges == 1) {
+ /*
+ * Two dots and an edge. This one's another
+ * quadratic equation.
+ *
+ * The point we want must lie on the perpendicular
+ * bisector of the two dots; that much is obvious.
+ * So we can construct a parametrisation of that
+ * bisecting line, giving linear formulae for x,y
+ * in terms of t. We can also express the distance
+ * from the edge as such a linear formula.
+ *
+ * Then we set that equal to the radius of the
+ * circle passing through the two points, which is
+ * a Pythagoras exercise; that gives rise to a
+ * quadratic in t, which we solve.
+ */
+ double xt[2], yt[2], rt[2]; /* a,b: {x,y,r}=at+b */
+ double q[3]; /* a,b,c: at^2+bt+c=0 */
+ double disc;
+ double halfsep;
+
+ /* Find parametric formulae for x,y. */
+ {
+ int x1 = dots[0]->x, x2 = dots[1]->x;
+ int y1 = dots[0]->y, y2 = dots[1]->y;
+ int dx = x2-x1, dy = y2-y1;
+ double d = sqrt((double)dx*dx + (double)dy*dy);
+
+ xt[1] = (x1+x2)/2.0;
+ yt[1] = (y1+y2)/2.0;
+ /* It's convenient if we have t at standard scale. */
+ xt[0] = -dy/d;
+ yt[0] = dx/d;
+
+ /* Also note down half the separation between
+ * the dots, for use in computing the circle radius. */
+ halfsep = 0.5*d;
+ }
+
+ /* Find a parametric formula for r. */
+ {
+ int x1 = edgedot1[0]->x, x2 = edgedot2[0]->x;
+ int y1 = edgedot1[0]->y, y2 = edgedot2[0]->y;
+ int dx = x2-x1, dy = y2-y1;
+ double d = sqrt((double)dx*dx + (double)dy*dy);
+ rt[0] = (xt[0]*dy - yt[0]*dx) / d;
+ rt[1] = ((xt[1]-x1)*dy - (yt[1]-y1)*dx) / d;
+ }
+
+ /* Construct the quadratic equation. */
+ q[0] = rt[0]*rt[0];
+ q[1] = 2*rt[0]*rt[1];
+ q[2] = rt[1]*rt[1];
+ q[0] -= 1;
+ q[2] -= halfsep*halfsep;
+
+ /* And solve it. */
+ disc = q[1]*q[1] - 4*q[0]*q[2];
+ if (disc >= 0) {
+ double t;
+
+ disc = sqrt(disc);
+
+ t = (-q[1] + disc) / (2*q[0]);
+ cx[cn] = xt[0]*t + xt[1];
+ cy[cn] = yt[0]*t + yt[1];
+ cn++;
+
+ t = (-q[1] - disc) / (2*q[0]);
+ cx[cn] = xt[0]*t + xt[1];
+ cy[cn] = yt[0]*t + yt[1];
+ cn++;
+ }
+ } else if (nedges == 0) {
+ /*
+ * Three dots. This is another linear matrix
+ * equation, this time with each row of the matrix
+ * representing the perpendicular bisector between
+ * two of the points. Of course we only need two
+ * such lines to find their intersection, so we
+ * need only solve a 2x2 matrix equation.
+ */
+
+ double matrix[4], vector[2], vector2[2];
+ int m;
+
+ for (m = 0; m < 2; m++) {
+ int x1 = dots[m]->x, x2 = dots[m+1]->x;
+ int y1 = dots[m]->y, y2 = dots[m+1]->y;
+ int dx = x2-x1, dy = y2-y1;
+
+ /*
+ * ((x,y) - (x1,y1)) . (dx,dy) = 1/2 |(dx,dy)|^2
+ *
+ * => 2x dx + 2y dy = dx^2+dy^2 + (2 x1 dx + 2 y1 dy)
+ */
+ matrix[2*m+0] = 2*dx;
+ matrix[2*m+1] = 2*dy;
+ vector[m] = ((double)dx*dx + (double)dy*dy +
+ 2.0*x1*dx + 2.0*y1*dy);
+ }
+
+ if (solve_2x2_matrix(matrix, vector, vector2)) {
+ cx[cn] = vector2[0];
+ cy[cn] = vector2[1];
+ cn++;
+ }
+ }
+
+ /*
+ * Now go through our candidate points and see if any
+ * of them are better than what we've got so far.
+ */
+ for (m = 0; m < cn; m++) {
+ double x = cx[m], y = cy[m];
+
+ /*
+ * First, disqualify the point if it's not inside
+ * the polygon, which we work out by counting the
+ * edges to the right of the point. (For
+ * tiebreaking purposes when edges start or end on
+ * our y-coordinate or go right through it, we
+ * consider our point to be offset by a small
+ * _positive_ epsilon in both the x- and
+ * y-direction.)
+ */
+ int e, in = 0;
+ for (e = 0; e < f->order; e++) {
+ int xs = f->edges[e]->dot1->x;
+ int xe = f->edges[e]->dot2->x;
+ int ys = f->edges[e]->dot1->y;
+ int ye = f->edges[e]->dot2->y;
+ if ((y >= ys && y < ye) || (y >= ye && y < ys)) {
+ /*
+ * The line goes past our y-position. Now we need
+ * to know if its x-coordinate when it does so is
+ * to our right.
+ *
+ * The x-coordinate in question is mathematically
+ * (y - ys) * (xe - xs) / (ye - ys), and we want
+ * to know whether (x - xs) >= that. Of course we
+ * avoid the division, so we can work in integers;
+ * to do this we must multiply both sides of the
+ * inequality by ye - ys, which means we must
+ * first check that's not negative.
+ */
+ int num = xe - xs, denom = ye - ys;
+ if (denom < 0) {
+ num = -num;
+ denom = -denom;
+ }
+ if ((x - xs) * denom >= (y - ys) * num)
+ in ^= 1;
+ }
+ }
+
+ if (in) {
+ double mindist = HUGE_VAL;
+ int e, d;
+
+ /*
+ * This point is inside the polygon, so now we check
+ * its minimum distance to every edge and corner.
+ * First the corners ...
+ */
+ for (d = 0; d < f->order; d++) {
+ int xp = f->dots[d]->x;
+ int yp = f->dots[d]->y;
+ double dx = x - xp, dy = y - yp;
+ double dist = dx*dx + dy*dy;
+ if (mindist > dist)
+ mindist = dist;
+ }
+
+ /*
+ * ... and now also check the perpendicular distance
+ * to every edge, if the perpendicular lies between
+ * the edge's endpoints.
+ */
+ for (e = 0; e < f->order; e++) {
+ int xs = f->edges[e]->dot1->x;
+ int xe = f->edges[e]->dot2->x;
+ int ys = f->edges[e]->dot1->y;
+ int ye = f->edges[e]->dot2->y;
+
+ /*
+ * If s and e are our endpoints, and p our
+ * candidate circle centre, the foot of a
+ * perpendicular from p to the line se lies
+ * between s and e if and only if (p-s).(e-s) lies
+ * strictly between 0 and (e-s).(e-s).
+ */
+ int edx = xe - xs, edy = ye - ys;
+ double pdx = x - xs, pdy = y - ys;
+ double pde = pdx * edx + pdy * edy;
+ long ede = (long)edx * edx + (long)edy * edy;
+ if (0 < pde && pde < ede) {
+ /*
+ * Yes, the nearest point on this edge is
+ * closer than either endpoint, so we must
+ * take it into account by measuring the
+ * perpendicular distance to the edge and
+ * checking its square against mindist.
+ */
+
+ double pdre = pdx * edy - pdy * edx;
+ double sqlen = pdre * pdre / ede;
+
+ if (mindist > sqlen)
+ mindist = sqlen;
+ }
+ }
+
+ /*
+ * Right. Now we know the biggest circle around this
+ * point, so we can check it against bestdist.
+ */
+ if (bestdist < mindist) {
+ bestdist = mindist;
+ xbest = x;
+ ybest = y;
+ }
+ }
+ }
+
+ if (k < f->order)
+ nedges--;
+ else
+ ndots--;
+ }
+ if (j < f->order)
+ nedges--;
+ else
+ ndots--;
+ }
+ if (i < f->order)
+ nedges--;
+ else
+ ndots--;
+ }
+
+ assert(bestdist > 0);
+
+ f->has_incentre = TRUE;
+ f->ix = xbest + 0.5; /* round to nearest */
+ f->iy = ybest + 0.5;
+}
+
/* ------ Generate various types of grid ------ */
/* General method is to generate faces, by calculating their dot coordinates.
freetree234(points);
assert(g->num_faces <= max_faces);
assert(g->num_dots <= max_dots);
- g->middle_face = g->faces + (height/2) * width + (width/2);
grid_make_consistent(g);
return g;
freetree234(points);
assert(g->num_faces <= max_faces);
assert(g->num_dots <= max_dots);
- g->middle_face = g->faces + (height/2) * width + (width/2);
grid_make_consistent(g);
return g;
}
}
- /* "+ width" takes us to the middle of the row, because each row has
- * (2*width) faces. */
- g->middle_face = g->faces + (height / 2) * 2 * width + width;
-
grid_make_consistent(g);
return g;
}
freetree234(points);
assert(g->num_faces <= max_faces);
assert(g->num_dots <= max_dots);
- g->middle_face = g->faces + (height/2) * width + (width/2);
grid_make_consistent(g);
return g;
freetree234(points);
assert(g->num_faces <= max_faces);
assert(g->num_dots <= max_dots);
- g->middle_face = g->faces + (height/2) * width + (width/2);
grid_make_consistent(g);
return g;
freetree234(points);
assert(g->num_faces <= max_faces);
assert(g->num_dots <= max_dots);
- g->middle_face = g->faces + (height/2) * width + (width/2);
grid_make_consistent(g);
return g;
freetree234(points);
assert(g->num_faces <= max_faces);
assert(g->num_dots <= max_dots);
- g->middle_face = g->faces + (height/2) * width + (width/2);
grid_make_consistent(g);
return g;
freetree234(points);
assert(g->num_faces <= max_faces);
assert(g->num_dots <= max_dots);
- g->middle_face = g->faces + 6 * ((height/2) * width + (width/2));
+
+ grid_make_consistent(g);
+ return g;
+}
+
+grid *grid_new_floret(int width, int height)
+{
+ int x, y;
+ /* Vectors for sides; weird numbers needed to keep puzzle aligned with window
+ * -py/px is close to tan(30 - atan(sqrt(3)/9))
+ * using py=26 makes everything lean to the left, rather than right
+ */
+ int px = 75, py = -26; /* |( 75, -26)| = 79.43 */
+ int qx = 4*px/5, qy = -py*2; /* |( 60, 52)| = 79.40 */
+ int rx = qx-px, ry = qy-py; /* |(-15, 78)| = 79.38 */
+
+ /* Upper bounds - don't have to be exact */
+ int max_faces = 6 * width * height;
+ int max_dots = 9 * (width + 1) * (height + 1);
+
+ tree234 *points;
+
+ grid *g = grid_new();
+ g->tilesize = 2 * px;
+ g->faces = snewn(max_faces, grid_face);
+ g->dots = snewn(max_dots, grid_dot);
+
+ points = newtree234(grid_point_cmp_fn);
+
+ /* generate pentagonal faces */
+ for (y = 0; y < height; y++) {
+ for (x = 0; x < width; x++) {
+ grid_dot *d;
+ /* face centre */
+ int cx = (6*px+3*qx)/2 * x;
+ int cy = (4*py-5*qy) * y;
+ if (x % 2)
+ cy -= (4*py-5*qy)/2;
+ else if (y && y == height-1)
+ continue; /* make better looking grids? try 3x3 for instance */
+
+ grid_face_add_new(g, 5);
+ d = grid_get_dot(g, points, cx , cy ); grid_face_set_dot(g, d, 0);
+ d = grid_get_dot(g, points, cx+2*rx , cy+2*ry ); grid_face_set_dot(g, d, 1);
+ d = grid_get_dot(g, points, cx+2*rx+qx, cy+2*ry+qy); grid_face_set_dot(g, d, 2);
+ d = grid_get_dot(g, points, cx+2*qx+rx, cy+2*qy+ry); grid_face_set_dot(g, d, 3);
+ d = grid_get_dot(g, points, cx+2*qx , cy+2*qy ); grid_face_set_dot(g, d, 4);
+
+ grid_face_add_new(g, 5);
+ d = grid_get_dot(g, points, cx , cy ); grid_face_set_dot(g, d, 0);
+ d = grid_get_dot(g, points, cx+2*qx , cy+2*qy ); grid_face_set_dot(g, d, 1);
+ d = grid_get_dot(g, points, cx+2*qx+px, cy+2*qy+py); grid_face_set_dot(g, d, 2);
+ d = grid_get_dot(g, points, cx+2*px+qx, cy+2*py+qy); grid_face_set_dot(g, d, 3);
+ d = grid_get_dot(g, points, cx+2*px , cy+2*py ); grid_face_set_dot(g, d, 4);
+
+ grid_face_add_new(g, 5);
+ d = grid_get_dot(g, points, cx , cy ); grid_face_set_dot(g, d, 0);
+ d = grid_get_dot(g, points, cx+2*px , cy+2*py ); grid_face_set_dot(g, d, 1);
+ d = grid_get_dot(g, points, cx+2*px-rx, cy+2*py-ry); grid_face_set_dot(g, d, 2);
+ d = grid_get_dot(g, points, cx-2*rx+px, cy-2*ry+py); grid_face_set_dot(g, d, 3);
+ d = grid_get_dot(g, points, cx-2*rx , cy-2*ry ); grid_face_set_dot(g, d, 4);
+
+ grid_face_add_new(g, 5);
+ d = grid_get_dot(g, points, cx , cy ); grid_face_set_dot(g, d, 0);
+ d = grid_get_dot(g, points, cx-2*rx , cy-2*ry ); grid_face_set_dot(g, d, 1);
+ d = grid_get_dot(g, points, cx-2*rx-qx, cy-2*ry-qy); grid_face_set_dot(g, d, 2);
+ d = grid_get_dot(g, points, cx-2*qx-rx, cy-2*qy-ry); grid_face_set_dot(g, d, 3);
+ d = grid_get_dot(g, points, cx-2*qx , cy-2*qy ); grid_face_set_dot(g, d, 4);
+
+ grid_face_add_new(g, 5);
+ d = grid_get_dot(g, points, cx , cy ); grid_face_set_dot(g, d, 0);
+ d = grid_get_dot(g, points, cx-2*qx , cy-2*qy ); grid_face_set_dot(g, d, 1);
+ d = grid_get_dot(g, points, cx-2*qx-px, cy-2*qy-py); grid_face_set_dot(g, d, 2);
+ d = grid_get_dot(g, points, cx-2*px-qx, cy-2*py-qy); grid_face_set_dot(g, d, 3);
+ d = grid_get_dot(g, points, cx-2*px , cy-2*py ); grid_face_set_dot(g, d, 4);
+
+ grid_face_add_new(g, 5);
+ d = grid_get_dot(g, points, cx , cy ); grid_face_set_dot(g, d, 0);
+ d = grid_get_dot(g, points, cx-2*px , cy-2*py ); grid_face_set_dot(g, d, 1);
+ d = grid_get_dot(g, points, cx-2*px+rx, cy-2*py+ry); grid_face_set_dot(g, d, 2);
+ d = grid_get_dot(g, points, cx+2*rx-px, cy+2*ry-py); grid_face_set_dot(g, d, 3);
+ d = grid_get_dot(g, points, cx+2*rx , cy+2*ry ); grid_face_set_dot(g, d, 4);
+ }
+ }
+
+ freetree234(points);
+ assert(g->num_faces <= max_faces);
+ assert(g->num_dots <= max_dots);
+
+ grid_make_consistent(g);
+ return g;
+}
+
+grid *grid_new_dodecagonal(int width, int height)
+{
+ int x, y;
+ /* Vector for side of triangle - ratio is close to sqrt(3) */
+ int a = 15;
+ int b = 26;
+
+ /* Upper bounds - don't have to be exact */
+ int max_faces = 3 * width * height;
+ int max_dots = 14 * width * height;
+
+ tree234 *points;
+
+ grid *g = grid_new();
+ g->tilesize = b;
+ g->faces = snewn(max_faces, grid_face);
+ g->dots = snewn(max_dots, grid_dot);
+
+ points = newtree234(grid_point_cmp_fn);
+
+ for (y = 0; y < height; y++) {
+ for (x = 0; x < width; x++) {
+ grid_dot *d;
+ /* centre of dodecagon */
+ int px = (4*a + 2*b) * x;
+ int py = (3*a + 2*b) * y;
+ if (y % 2)
+ px += 2*a + b;
+
+ /* dodecagon */
+ grid_face_add_new(g, 12);
+ d = grid_get_dot(g, points, px + ( a ), py - (2*a + b)); grid_face_set_dot(g, d, 0);
+ d = grid_get_dot(g, points, px + ( a + b), py - ( a + b)); grid_face_set_dot(g, d, 1);
+ d = grid_get_dot(g, points, px + (2*a + b), py - ( a )); grid_face_set_dot(g, d, 2);
+ d = grid_get_dot(g, points, px + (2*a + b), py + ( a )); grid_face_set_dot(g, d, 3);
+ d = grid_get_dot(g, points, px + ( a + b), py + ( a + b)); grid_face_set_dot(g, d, 4);
+ d = grid_get_dot(g, points, px + ( a ), py + (2*a + b)); grid_face_set_dot(g, d, 5);
+ d = grid_get_dot(g, points, px - ( a ), py + (2*a + b)); grid_face_set_dot(g, d, 6);
+ d = grid_get_dot(g, points, px - ( a + b), py + ( a + b)); grid_face_set_dot(g, d, 7);
+ d = grid_get_dot(g, points, px - (2*a + b), py + ( a )); grid_face_set_dot(g, d, 8);
+ d = grid_get_dot(g, points, px - (2*a + b), py - ( a )); grid_face_set_dot(g, d, 9);
+ d = grid_get_dot(g, points, px - ( a + b), py - ( a + b)); grid_face_set_dot(g, d, 10);
+ d = grid_get_dot(g, points, px - ( a ), py - (2*a + b)); grid_face_set_dot(g, d, 11);
+
+ /* triangle below dodecagon */
+ if ((y < height - 1 && (x < width - 1 || !(y % 2)) && (x > 0 || (y % 2)))) {
+ grid_face_add_new(g, 3);
+ d = grid_get_dot(g, points, px + a, py + (2*a + b)); grid_face_set_dot(g, d, 0);
+ d = grid_get_dot(g, points, px , py + (2*a + 2*b)); grid_face_set_dot(g, d, 1);
+ d = grid_get_dot(g, points, px - a, py + (2*a + b)); grid_face_set_dot(g, d, 2);
+ }
+
+ /* triangle above dodecagon */
+ if ((y && (x < width - 1 || !(y % 2)) && (x > 0 || (y % 2)))) {
+ grid_face_add_new(g, 3);
+ d = grid_get_dot(g, points, px - a, py - (2*a + b)); grid_face_set_dot(g, d, 0);
+ d = grid_get_dot(g, points, px , py - (2*a + 2*b)); grid_face_set_dot(g, d, 1);
+ d = grid_get_dot(g, points, px + a, py - (2*a + b)); grid_face_set_dot(g, d, 2);
+ }
+ }
+ }
+
+ freetree234(points);
+ assert(g->num_faces <= max_faces);
+ assert(g->num_dots <= max_dots);
+
+ grid_make_consistent(g);
+ return g;
+}
+
+grid *grid_new_greatdodecagonal(int width, int height)
+{
+ int x, y;
+ /* Vector for side of triangle - ratio is close to sqrt(3) */
+ int a = 15;
+ int b = 26;
+
+ /* Upper bounds - don't have to be exact */
+ int max_faces = 30 * width * height;
+ int max_dots = 200 * width * height;
+
+ tree234 *points;
+
+ grid *g = grid_new();
+ g->tilesize = b;
+ g->faces = snewn(max_faces, grid_face);
+ g->dots = snewn(max_dots, grid_dot);
+
+ points = newtree234(grid_point_cmp_fn);
+
+ for (y = 0; y < height; y++) {
+ for (x = 0; x < width; x++) {
+ grid_dot *d;
+ /* centre of dodecagon */
+ int px = (6*a + 2*b) * x;
+ int py = (3*a + 3*b) * y;
+ if (y % 2)
+ px += 3*a + b;
+
+ /* dodecagon */
+ grid_face_add_new(g, 12);
+ d = grid_get_dot(g, points, px + ( a ), py - (2*a + b)); grid_face_set_dot(g, d, 0);
+ d = grid_get_dot(g, points, px + ( a + b), py - ( a + b)); grid_face_set_dot(g, d, 1);
+ d = grid_get_dot(g, points, px + (2*a + b), py - ( a )); grid_face_set_dot(g, d, 2);
+ d = grid_get_dot(g, points, px + (2*a + b), py + ( a )); grid_face_set_dot(g, d, 3);
+ d = grid_get_dot(g, points, px + ( a + b), py + ( a + b)); grid_face_set_dot(g, d, 4);
+ d = grid_get_dot(g, points, px + ( a ), py + (2*a + b)); grid_face_set_dot(g, d, 5);
+ d = grid_get_dot(g, points, px - ( a ), py + (2*a + b)); grid_face_set_dot(g, d, 6);
+ d = grid_get_dot(g, points, px - ( a + b), py + ( a + b)); grid_face_set_dot(g, d, 7);
+ d = grid_get_dot(g, points, px - (2*a + b), py + ( a )); grid_face_set_dot(g, d, 8);
+ d = grid_get_dot(g, points, px - (2*a + b), py - ( a )); grid_face_set_dot(g, d, 9);
+ d = grid_get_dot(g, points, px - ( a + b), py - ( a + b)); grid_face_set_dot(g, d, 10);
+ d = grid_get_dot(g, points, px - ( a ), py - (2*a + b)); grid_face_set_dot(g, d, 11);
+
+ /* hexagon below dodecagon */
+ if (y < height - 1 && (x < width - 1 || !(y % 2)) && (x > 0 || (y % 2))) {
+ grid_face_add_new(g, 6);
+ d = grid_get_dot(g, points, px + a, py + (2*a + b)); grid_face_set_dot(g, d, 0);
+ d = grid_get_dot(g, points, px + 2*a, py + (2*a + 2*b)); grid_face_set_dot(g, d, 1);
+ d = grid_get_dot(g, points, px + a, py + (2*a + 3*b)); grid_face_set_dot(g, d, 2);
+ d = grid_get_dot(g, points, px - a, py + (2*a + 3*b)); grid_face_set_dot(g, d, 3);
+ d = grid_get_dot(g, points, px - 2*a, py + (2*a + 2*b)); grid_face_set_dot(g, d, 4);
+ d = grid_get_dot(g, points, px - a, py + (2*a + b)); grid_face_set_dot(g, d, 5);
+ }
+
+ /* hexagon above dodecagon */
+ if (y && (x < width - 1 || !(y % 2)) && (x > 0 || (y % 2))) {
+ grid_face_add_new(g, 6);
+ d = grid_get_dot(g, points, px - a, py - (2*a + b)); grid_face_set_dot(g, d, 0);
+ d = grid_get_dot(g, points, px - 2*a, py - (2*a + 2*b)); grid_face_set_dot(g, d, 1);
+ d = grid_get_dot(g, points, px - a, py - (2*a + 3*b)); grid_face_set_dot(g, d, 2);
+ d = grid_get_dot(g, points, px + a, py - (2*a + 3*b)); grid_face_set_dot(g, d, 3);
+ d = grid_get_dot(g, points, px + 2*a, py - (2*a + 2*b)); grid_face_set_dot(g, d, 4);
+ d = grid_get_dot(g, points, px + a, py - (2*a + b)); grid_face_set_dot(g, d, 5);
+ }
+
+ /* square on right of dodecagon */
+ if (x < width - 1) {
+ grid_face_add_new(g, 4);
+ d = grid_get_dot(g, points, px + 2*a + b, py - a); grid_face_set_dot(g, d, 0);
+ d = grid_get_dot(g, points, px + 4*a + b, py - a); grid_face_set_dot(g, d, 1);
+ d = grid_get_dot(g, points, px + 4*a + b, py + a); grid_face_set_dot(g, d, 2);
+ d = grid_get_dot(g, points, px + 2*a + b, py + a); grid_face_set_dot(g, d, 3);
+ }
+
+ /* square on top right of dodecagon */
+ if (y && (x < width - 1 || !(y % 2))) {
+ grid_face_add_new(g, 4);
+ d = grid_get_dot(g, points, px + ( a ), py - (2*a + b)); grid_face_set_dot(g, d, 0);
+ d = grid_get_dot(g, points, px + (2*a ), py - (2*a + 2*b)); grid_face_set_dot(g, d, 1);
+ d = grid_get_dot(g, points, px + (2*a + b), py - ( a + 2*b)); grid_face_set_dot(g, d, 2);
+ d = grid_get_dot(g, points, px + ( a + b), py - ( a + b)); grid_face_set_dot(g, d, 3);
+ }
+
+ /* square on top left of dodecagon */
+ if (y && (x || (y % 2))) {
+ grid_face_add_new(g, 4);
+ d = grid_get_dot(g, points, px - ( a + b), py - ( a + b)); grid_face_set_dot(g, d, 0);
+ d = grid_get_dot(g, points, px - (2*a + b), py - ( a + 2*b)); grid_face_set_dot(g, d, 1);
+ d = grid_get_dot(g, points, px - (2*a ), py - (2*a + 2*b)); grid_face_set_dot(g, d, 2);
+ d = grid_get_dot(g, points, px - ( a ), py - (2*a + b)); grid_face_set_dot(g, d, 3);
+ }
+ }
+ }
+
+ freetree234(points);
+ assert(g->num_faces <= max_faces);
+ assert(g->num_dots <= max_dots);
grid_make_consistent(g);
return g;
~mdw / sgt / puzzles / blobdiff