From 8646c6b2ddaf32ac7342cc8b61e559d46885af4f Mon Sep 17 00:00:00 2001 From: root <> Date: Sun, 8 Feb 2009 16:49:14 +0000 Subject: *** empty log message *** --- src/contour.c | 383 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 383 insertions(+) create mode 100644 src/contour.c (limited to 'src/contour.c') diff --git a/src/contour.c b/src/contour.c new file mode 100644 index 0000000..8022adb --- /dev/null +++ b/src/contour.c @@ -0,0 +1,383 @@ +#include "project.h" + +#ifndef PI +#define PI 3.141592654 +#endif + + +typedef struct { + int x, y; +} Point; + +typedef struct contour_struct { + float *data; + float **dcol; + unsigned char *flags; + unsigned char **fcol; + int x, y; + int x1, y1; + Jwgline current_line; + int current_line_size; + int lp; + Handle h; +} *Contour; + + +#define INRANGE(l1,l2,c) ( ((c)<(l2)) && ((c)>=(l1)) ) +#define RANGEX(c,x) ( ((x)>0) ? ( ((x)<((c)->x1)) ? (x):((c)->x1)) : 1) +#define RANGEY(c,y) (((y)>0) ? ( ((y)<((c)->y1)) ? (y):((c)->y1)) : 1) + +#define FLAGS(c,p) ((c->fcol[p.y])[p.x]) +#define DATA(c,p) ((c->dcol[p.y])[p.x]) + +#define DEF_LINE_SIZE 8192 + + + + + +void +cont_reset_point(Contour c) +{ + c->current_line.npts = 0; +} + +void +cont_add_point(Contour c, Jwgpos pos) +{ + Jwgpos *ptr; + int n; + + if (c->current_line.npts == c->current_line_size) { + /* uh-oh we're about to over flow */ + + ptr = c->current_line.data; + + + n = c->current_line_size; + printf("JWG: CONTOUR: Increasing maximum line size from %d to %d\n", + n, n * 2); + + c->current_line_size = n * 2; + + c->current_line.data = malloc(sizeof(Jwgpos) * n * 2); + bcopy(ptr, c->current_line.data, sizeof(Jwgpos) * n); + + free(ptr); + + } + jwg_xform(c->h, &pos); + + c->current_line.data[c->current_line.npts] = pos; + c->current_line.npts++; + +} + +unsigned char +find_direction(Point p1, Point p2) +{ + p2.x -= p1.x; + p2.y -= p1.y; + + switch (p2.y) { + case -1: + switch (p2.x) { + case -1: + return (1); + case 0: + return (2); + case 1: + return (4); + } + case 0: + switch (p2.x) { + case -1: + return (128); + case 1: + return (8); + } + case 1: + switch (p2.x) { + case -1: + return (64); + case 0: + return (32); + case 1: + return (16); + } + } + + + fprintf(stderr, "Contouring cockup type 3 ... (bye bye)\n"); + abort(); + +} + +Point +add_direction(Point p, int d) +{ + Point dirs[8] = {{-1, -1}, {0, -1}, {1, -1}, {1, 0}, {1, 1}, {0, 1}, {-1, 1}, {-1, 0}}; + + p.x = dirs[d].x; + p.y = dirs[d].y; + + return (p); +} + + +int +points_marked(Contour c, Point p1, Point p2) +{ + return (c->fcol[p1.y])[p1.x] & find_direction(p1, p2); +} + +void +mark_points(Contour c, Point p1, Point p2) +{ + (c->fcol[p1.y])[p1.x] |= find_direction(p1, p2); + (c->fcol[p2.y])[p2.x] |= find_direction(p2, p1); +} + + + +trace(Contour c, Point p1, Point p2, Point pl, float level) +{ + Point pn; + float v1, v2, vn; + Jwgpos pos; + float r; + +printf("."); +fflush(stdout); + + if (!c->lp) + { + pos.x = 1; + pos.y = 1; + + cont_add_point(c, pos); + } + + v1 = DATA(c, p1); + v2 = DATA(c, p2); + + while (1) { + + /* + * printf("(%d,%d,%f) (%d,%d,%f) (%d,%d,%f)\n", p1.x, p1.y, DATA(c, p1), + * p2.x, p2.y, DATA(c, p2), pl.x, pl.y, DATA(c, pl)); + */ + + r = (v2 - v1); + if (r == 0.0) { + r = 0.5; + } else { + r = (level - v1) / r; + } + + pos.x = r * (float) RANGEX(c, p2.x); + pos.y = r * (float) RANGEY(c, p2.y); + + r = 1.0 - r; + + pos.x += r * (float) RANGEX(c, p1.x); + pos.y += r * (float) RANGEY(c, p1.y); + + cont_add_point(c, pos); + + + + + + + + if (DATA(c, p1) > level) { + printf("Contouring cockup type 1 -- bye bye\n"); + abort(); + } + if (DATA(c, p2) <= level) { + printf("Contouring cockup type 2 -- bye bye\n"); + abort(); + } + if (points_marked(c, p1, p2)) { + if (c->lp) { + jwg_write_line(c->h, &c->current_line); + c->current_line.npts = 0; + } else { + pos.x = 1; + pos.y = 1; + + cont_add_point(c, pos); + } + + + return; + } + mark_points(c, p1, p2); + + pn.x = p1.x + p2.x - pl.x; + pn.y = p1.y + p2.y - pl.y; + + vn = DATA(c, pn); + + if (vn <= level) { + pl = p1; + p1 = pn; + v1 = vn; + } else { + pl = p2; + p2 = pn; + v2 = vn; + } + + + + + } + +} + + + + + +possible(Contour c, int x1, int y1, int x2, int y2, int xl, int yl, float level) +{ + Point p1, p2, pl; + + p1.x = x1; + p1.y = y1; + p2.x = x2; + p2.y = y2; + pl.x = xl; + pl.y = yl; + + if (points_marked(c, p1, p2)) + return; + + trace(c, p1, p2, pl, level); + + +} + + +contour(Contour c, float level) +{ + int x, y; + float **c0, **c1; + float *r00, *r10, *r11, *r01; + + c1 = c0 = c->dcol; + c1--; + + for (y = 0; y < c->y; ++y) { + r00 = r01 = *(c0); + r01--; + if (y) { + r10 = r11 = *(c1); + r11--; + } + for (x = 0; x < c->x; ++x) { + + if (x) { + if (INRANGE(*r00, *r01, level)) + possible(c, x, y, x - 1, y, x, y - 1, level); + if (INRANGE(*r01, *r00, level)) + possible(c, x - 1, y, x, y, x - 1, y + 1, level); + } + if (y) { + if (INRANGE(*r00, *r10, level)) + possible(c, x, y - 1, x, y, x - 1, y, level); + if (INRANGE(*r10, *r00, level)) + possible(c, x, y, x, y - 1, x + 1, y - 1, level); + } + r00++; + r10++; + r11++; + r01++; + } + c0++; + c1++; + } + +} + + + +void +do_contour(Handle h, float level, int side, int lp) +{ + Contour c = (Contour) malloc(sizeof(struct contour_struct)); + int i; + float hug = 1.0e10 * (float) side; + + + c->data = h->data->data; + c->x1 = c->x = h->data->w; + c->y1 = c->y = h->data->h; + + + c->dcol = (float **) malloc(sizeof(float *) * c->y); + + + + c->dcol[0] = c->data; + for (i = 1; i < c->y; ++i) + c->dcol[i] = c->dcol[i - 1] + c->x; + + c->flags = (unsigned char *) malloc(c->x * c->y); + bzero(c->flags, c->x * c->y); + + c->fcol = (unsigned char **) malloc(c->y * sizeof(unsigned char *)); + + c->fcol[0] = c->flags; + for (i = 1; i < c->y; ++i) + c->fcol[i] = c->fcol[i - 1] + c->x; + + c->x1--; + c->y1--; + + for (i = 0; i < c->y; ++i) { + (c->dcol[i])[0] = hug; + (c->dcol[i])[c->x1] = hug; + } + for (i = 0; i < c->x; ++i) { + (c->dcol[0])[i] = hug; + (c->dcol[c->y1])[i] = hug; + } + + c->x1--; + c->y1--; + + /* + * { Point p; + * + * for (p.y=0;p.yy;++p.y) { for (p.x=0;p.xx;++p.x) { + * + * printf("%.1f ",DATA(c,p)); + * + * } printf("\n"); + * + * } + * } */ + + c->current_line.npts = 0; + c->current_line_size = DEF_LINE_SIZE; + c->current_line.data = (Jwgpos *) malloc(sizeof(Jwgpos) * DEF_LINE_SIZE); + + c->h = h; + c->lp = lp; + + contour(c, level); + printf("\n"); + + + if (!c->lp) + jwg_write_polygon(c->h, &c->current_line); + + free(c->flags); + free(c->fcol); + free(c->dcol); + free(c->current_line.data); + +} -- cgit v1.2.3