1 /* graph.c 1.10 84/03/03 2 * 3 * This file contains the functions for producing the graphics 4 * images in the varian/versatec drivers for ditroff. 5 */ 6 7 8 #include <stdio.h> 9 #include <ctype.h> 10 #include <math.h> 11 12 13 #define TRUE 1 14 #define FALSE 0 15 /* imports from dver.c */ 16 #define hmot(n) hpos += n; 17 #define vmot(n) vgoto(vpos + n); 18 19 extern int hpos; 20 extern int vpos; 21 extern int output; 22 extern vgoto(); 23 extern point(); 24 25 #define MAXPOINTS 200 /* number of points legal for a curve */ 26 27 #define SOLID -1 /* line styles: these are used as bit masks to */ 28 #define DOTTED 004 /* create the right style lines. */ 29 #define DASHED 020 30 #define DOTDASHED 024 31 #define LONGDASHED 074 32 /* constants... */ 33 #define pi 3.14159265358979324 34 #define log2_10 3.3219280948873623 35 36 37 int linethickness = 3; /* number of pixels wide to make lines */ 38 int linmod = SOLID; /* type of line (SOLID, DOTTED, DASHED...) */ 39 40 41 42 /*---------------------------------------------------------------------------- 43 * Routine: drawline (horizontal_motion, vertical_motion) 44 * 45 * Results: Draws a line of "linethickness" width and "linmod" style 46 * from current (hpos, vpos) to (hpos + dh, vpos + dv). 47 * 48 * Side Efct: Resulting position is at end of line (hpos + dh, vpos + dv) 49 *----------------------------------------------------------------------------*/ 50 51 drawline(dh, dv) 52 register int dh; 53 register int dv; 54 { 55 if (output) HGtline (hpos, vpos, hpos + dh, vpos + dv); 56 hmot (dh); /* new position is at */ 57 vmot (dv); /* the end of the line */ 58 } 59 60 61 /*---------------------------------------------------------------------------- 62 * Routine: drawcirc (diameter) 63 * 64 * Results: Draws a circle with leftmost point at current (hpos, vpos) 65 * with the given diameter d. 66 * 67 * Side Efct: Resulting position is at (hpos + diameter, vpos) 68 *----------------------------------------------------------------------------*/ 69 70 drawcirc(d) 71 register int d; 72 { /* 0.0 is the angle to sweep the arc: = full circle */ 73 if (output) HGArc (hpos + d/2, vpos, hpos, vpos, 0.0); 74 hmot (d); /* new postion is the right of the circle */ 75 } 76 77 78 /*---------------------------------------------------------------------------- 79 * Routine: drawellip (horizontal_diameter, vertical_diameter) 80 * 81 * Results: Draws regular ellipses given the major "diameters." It does 82 * so using a modified circle algorithm (see RoundEnd) that 83 * increments x and y proportionally to their axes. 84 * 85 * Side Efct: Resulting position is at (hpos + hd, vpos). 86 *----------------------------------------------------------------------------*/ 87 88 drawellip(hd, vd) 89 int hd; 90 int vd; 91 { 92 double xs, ys, xepsilon, yepsilon; 93 register int thick; 94 register int basex; 95 register int basey; 96 register int x; 97 register int y; 98 99 100 basex = hpos + (hd >> 1); /* bases == coordinates of center */ 101 hmot(hd); /* horizontal motion (hpos should */ 102 basey = vpos; /* NOT be used after this) */ 103 /* hd and vd are radii, not diameters */ 104 if ((hd = hd >> 1) < 1) hd++; /* neither diameter can be zero. */ 105 if ((vd = vd >> 1) < 1) vd++; /* hd changed!! no hmot after this */ 106 ys = (double) vd; /* start at top of the ellipse */ 107 xs = 0.0; /* (y = 1/2 diameter, x = 0) */ 108 109 if ((thick = vd) > hd) thick = hd; 110 xepsilon = (double) thick / (double) (vd * vd); 111 yepsilon = (double) thick / (double) (hd * hd); 112 113 /* Calculate trajectory of the ellipse for 1/4 */ 114 /* the circumference (while ys is non-negative) */ 115 /* and mirror to get the other three quadrants. */ 116 117 thick = linethickness / 2; 118 if (thick) { /* more than one pixel thick */ 119 RoundEnd(basex, ((int)(ys + 0.5)) + basey, thick, 0); 120 RoundEnd(basex, basey - ((int)(ys + 0.5)), thick, 0); 121 while (ys >= 0) { 122 xs += xepsilon * ys; /* generate circumference */ 123 ys -= yepsilon * xs; 124 x = (int)(xs + 0.5); 125 y = (int)(ys + 0.5); 126 RoundEnd(x + basex, y + basey, thick, 0); 127 RoundEnd(x + basex, basey - y, thick, 0); 128 RoundEnd(basex - x, y + basey, thick, 0); 129 RoundEnd(basex - x, basey - y, thick, 0); 130 } 131 } else { /* do the perimeter only (no fill) */ 132 point(basex, ((int)(ys + 0.5)) + basey); 133 point(basex, basey - ((int)(ys + 0.5))); 134 while (ys >= 0) { 135 xs += xepsilon * ys; /* generate circumference */ 136 ys -= yepsilon * xs; 137 x = (int)(xs + 0.5); 138 y = (int)(ys + 0.5); 139 point(x + basex, y + basey); 140 point(x + basex, basey - y); 141 point(basex - x, y + basey); 142 point(basex - x, basey - y); 143 } 144 } 145 } 146 147 148 /*---------------------------------------------------------------------------- 149 * Routine: drawarc (xcenter, ycenter, xpoint, ypoint) 150 * 151 * Results: Draws an arc starting at current (hpos, vpos). Center is 152 * at (hpos + cdh, vpos + cdv) and the terminating point is 153 * at <center> + (pdh, pdv). The angle between the lines 154 * formed by the starting, ending, and center points is figured 155 * first, then the points and angle are sent to HGArc for the 156 * drawing. 157 * 158 * Side Efct: Resulting position is at the last point of the arc. 159 *----------------------------------------------------------------------------*/ 160 161 drawarc (cdh, cdv, pdh, pdv) 162 register int cdh; 163 register int cdv; 164 register int pdh; 165 register int pdv; 166 { 167 register double angle; 168 /* figure angle from the three points...*/ 169 /* and convert (and round) to degrees */ 170 angle = (atan2((double) pdh, (double) pdv) 171 - atan2 ((double) -cdh, (double) -cdv)) * 180.0 / pi; 172 /* "normalize" and round */ 173 angle += (angle < 0.0) ? 360.5 : 0.5; 174 175 if (output) HGArc(hpos + cdh, vpos + cdv, hpos, vpos, (int) angle); 176 hmot(cdh + pdh); 177 vmot(cdv + pdv); 178 } 179 180 181 /*---------------------------------------------------------------------------- 182 * Routine: drawwig (character_buffer, file_pointer, type_flag) 183 * 184 * Results: Given the starting position, the motion list in buf, and any 185 * extra characters from fp (terminated by a \n), drawwig sets 186 * up a point list to make a spline from. If "pic" is set picurve 187 * is called to draw the curve in pic style; else it calls HGCurve 188 * for the gremlin-style curve. 189 * 190 * Side Efct: Resulting position is reached from adding successive motions 191 * to the current position. 192 *----------------------------------------------------------------------------*/ 193 194 drawwig (buf, fp, pic) 195 char *buf; 196 FILE *fp; 197 int pic; 198 { 199 register int len = strlen(buf); /* length of the string in "buf" */ 200 register int npts = 2; /* point list index */ 201 register char *ptr = buf; /* "walking" pointer into buf */ 202 int x[MAXPOINTS], y[MAXPOINTS]; /* point list */ 203 204 while (*ptr == ' ') ptr++; /* skip any leading spaces */ 205 x[1] = hpos; /* the curve starts at the */ 206 y[1] = vpos; /* current position */ 207 208 while (*ptr != '\n') { /* curve commands end with a '\n' */ 209 hmot(atoi(ptr)); /* convert motion to curve points */ 210 x[npts] = hpos; /* and remember them */ 211 while (isdigit(*++ptr)); /* skip number*/ 212 while (*++ptr == ' '); /* skip spaces 'tween numbers */ 213 vmot(atoi(ptr)); 214 y[npts] = vpos; 215 while (isdigit(*++ptr)); 216 while (*ptr == ' ') ptr++; 217 /* if the amount we read wasn't the */ 218 /* whole thing, read some more in */ 219 if (len - (ptr - buf) < 15 && *(buf + len - 1) != '\n') { 220 char *cop = buf; 221 222 while (*cop++ = *ptr++); /* copy what's left to the beginning */ 223 fgets ((cop - 1), len - (cop - buf), fp); 224 ptr = buf; 225 } 226 if (npts < MAXPOINTS - 1) /* if too many points, forget some */ 227 npts++; 228 } 229 npts--; /* npts must point to the last coordinate in x and y */ 230 /* now, actually DO the curve */ 231 if (output) { 232 if (pic) 233 picurve(x, y, npts); 234 else 235 HGCurve(x, y, npts); 236 } 237 } 238 239 240 /*----------------------------------------------------------------------------* 241 | Routine: drawthick (thickness) 242 | 243 | Results: sets the variable "linethickness" to the given size. 244 | NO motion is involved. 245 *----------------------------------------------------------------------------*/ 246 247 drawthick(s) 248 int s; 249 { 250 linethickness = s; 251 } 252 253 254 /*----------------------------------------------------------------------------* 255 | Routine: drawstyle (style_bit_map) 256 | 257 | Results: sets the variable "linmod" to the given bit map. 258 | NO motion is involved. 259 *----------------------------------------------------------------------------*/ 260 261 drawstyle(s) 262 int s; 263 { 264 linmod = s; 265 } 266 267 268 /*---------------------------------------------------------------------------- 269 * Routine: picurve (xpoints, ypoints, num_of_points) 270 * 271 * Results: Draws a curve delimited by (not through) the line segments 272 * traced by (xpoints, ypoints) point list. This is the "Pic" 273 * style curve. 274 *----------------------------------------------------------------------------*/ 275 276 picurve (x, y, npts) 277 int x[MAXPOINTS]; 278 int y[MAXPOINTS]; 279 int npts; 280 { 281 register int i; /* line segment traverser */ 282 register int nseg; /* effective resolution for each curve */ 283 register float w; /* position factor */ 284 register int xp; /* current point (and intermediary) */ 285 register int yp; 286 int pxp, pyp; /* previous point (to make lines from) */ 287 float t1, t2, t3; /* calculation temps */ 288 int j; /* inner curve segment traverser */ 289 290 291 if (x[1] == x[npts] && y[1] == y[npts]) { 292 x[0] = x[npts - 1]; /* if the lines' ends meet, make */ 293 y[0] = y[npts - 1]; /* sure the curve meets */ 294 x[npts + 1] = x[2]; 295 y[npts + 1] = y[2]; 296 } else { /* otherwise, make the ends of the */ 297 x[0] = x[1]; /* curve touch the ending points of */ 298 y[0] = y[1]; /* the line segments */ 299 x[npts + 1] = x[npts]; 300 y[npts + 1] = y[npts]; 301 } 302 303 pxp = (x[0] + x[1]) / 2; /* make the last point pointers */ 304 pyp = (y[0] + y[1]) / 2; /* point to the start of the 1st line */ 305 306 for (i = 0; i < npts; i++) { /* traverse the line segments */ 307 xp = x[i] - x[i+1]; 308 yp = y[i] - y[i+1]; 309 nseg = (int) sqrt((double)(xp * xp + yp * yp)); 310 xp = x[i+1] - x[i+2]; 311 yp = y[i+1] - y[i+2]; /* "nseg" is the number of line */ 312 /* segments that will be drawn for */ 313 /* each curve segment. ">> 3" is */ 314 /* dropping the resolution ( == / 8) */ 315 nseg = (nseg + (int) sqrt((double)(xp * xp + yp * yp))) >> 3; 316 317 for (j = 1; j < nseg; j++) { 318 w = (float) j / (float) nseg; 319 t1 = 0.5 * w * w; 320 w -= 0.5; 321 t2 = 0.75 - w * w ; 322 w -= 0.5; 323 t3 = 0.5 * w * w; 324 xp = t1 * x[i+2] + t2 * x[i+1] + t3 * x[i] + 0.5; 325 yp = t1 * y[i+2] + t2 * y[i+1] + t3 * y[i] + 0.5; 326 327 HGtline(pxp, pyp, xp, yp); 328 pxp = xp; 329 pyp = yp; 330 } 331 } 332 } 333 334 335 /*---------------------------------------------------------------------------- 336 * Routine: line (xstart, ystart, xend, yend) 337 * 338 * Results: Draws a one-pixel wide line from (x0, y0) to (x1, y1) 339 * using point(x,y) to place the dots. Line style 340 * is done in this routine by masking off unwanted dots. 341 *----------------------------------------------------------------------------*/ 342 343 line(x0, y0, x1, y1) 344 register int x0; 345 register int y0; 346 int x1; 347 int y1; 348 { 349 register int dx; 350 register int dy; 351 register int res1; 352 int res2; 353 int xinc; 354 int yinc; 355 int slope; 356 357 xinc = 1; 358 yinc = 1; 359 if ((dx = x1-x0) < 0) { 360 xinc = -1; 361 dx = -dx; 362 } 363 if ((dy = y1-y0) < 0) { 364 yinc = -1; 365 dy = -dy; 366 } 367 slope = xinc*yinc; 368 res1 = 0; 369 res2 = 0; 370 if (dx >= dy) 371 while (x0 != x1) { 372 if((x0+slope*y0)&linmod) point(x0, y0); 373 if (res1 > res2) { 374 res2 += dx - res1; 375 res1 = 0; 376 y0 += yinc; 377 } 378 res1 += dy; 379 x0 += xinc; 380 } 381 else 382 while (y0 != y1) { 383 if((x0+slope*y0)&linmod) point(x0, y0); 384 if (res1 > res2) { 385 res2 += dy - res1; 386 res1 = 0; 387 x0 += xinc; 388 } 389 res1 += dx; 390 y0 += yinc; 391 } 392 if((x1+slope*y1)&linmod) point(x1, y1); 393 } 394 395 396 /*---------------------------------------------------------------------------- 397 * Routine: HGArc (xcenter, ycenter, xstart, ystart, angle) 398 * 399 * Results: This routine plots an arc centered about (cx, cy) counter 400 * clockwise starting from the point (px, py) through 'angle' 401 * degrees. If angle is 0, a full circle is drawn. 402 * It does so by calling RoundEnd (fat point maker) for points 403 * along the circle with density depending on the circle's size. 404 *----------------------------------------------------------------------------*/ 405 406 HGArc(cx,cy,px,py,angle) 407 register int cx; 408 register int cy; 409 int px; 410 int py; 411 int angle; 412 { 413 double xs, ys, resolution, epsilon, fullcircle; 414 int thick = linethickness / 2; 415 register int nx; 416 register int ny; 417 register int extent; 418 419 xs = px - cx; 420 ys = py - cy; 421 422 /* calculate drawing parameters */ 423 424 resolution = log10(sqrt( xs * xs + ys * ys)) * log2_10; 425 resolution = ceil(resolution); 426 resolution = pow(2.0, resolution); 427 epsilon = 1.0 / resolution; 428 fullcircle = 2 * pi * resolution; 429 fullcircle = ceil(fullcircle); 430 431 if (angle == 0) 432 extent = fullcircle; 433 else 434 extent = angle * fullcircle / 360.0; 435 436 while (extent-- > 0) { 437 xs += epsilon * ys; 438 nx = cx + (int) (xs + 0.5); 439 ys -= epsilon * xs; 440 ny = cy + (int) (ys + 0.5); 441 RoundEnd(nx, ny, thick, FALSE); 442 } /* end for */ 443 } /* end HGArc */ 444 445 446 /*---------------------------------------------------------------------------- 447 * Routine: RoundEnd (x, y, radius, filled_flag) 448 * 449 * Results: Plots a filled (if requested) circle of the specified radius 450 * centered about (x, y). 451 *----------------------------------------------------------------------------*/ 452 453 RoundEnd(x, y, radius, filled) 454 register int x; 455 register int y; 456 int radius; 457 int filled; 458 { 459 double xs, ys; /* floating point distance form center of circle */ 460 double epsilon; /* "resolution" of the step around circle */ 461 register int cx; /* center of circle */ 462 register int cy; 463 register int nx; /* integer distance from center */ 464 register int ny; 465 466 467 if (radius < 1) { /* too small to notice */ 468 point(x, y); 469 return; 470 } 471 472 xs = 0; 473 ys = radius; 474 epsilon = 1.0 / radius; 475 476 /* Calculate the trajectory of the circle for 1/4 the circumference 477 * and mirror appropriately to get the other three quadrants. 478 */ 479 480 nx = x; /* must start out the x and y for first */ 481 ny = y + radius; /* painting going on in while loop */ 482 483 while (ny >= 0) 484 { 485 if (filled) { /* fill from center */ 486 cx = x; 487 cy = y; 488 } else { /* fill from perimeter only (no fill) */ 489 cx = nx; 490 cy = ny; 491 } 492 while (cx <= nx) { 493 while (cy <= ny) { 494 point(cx, cy); 495 point(cx, 2*y-cy); 496 point(2*x-cx, cy); 497 point(2*x-cx, 2*y-cy); 498 cy++; 499 } /* end for k */ 500 cx++; 501 } /* end for j */; 502 /* generate circumference */ 503 xs += epsilon * ys; 504 nx = x + (int) (xs + 0.5); 505 ys -= epsilon * xs; 506 ny = y + (int) (ys + 0.5); 507 } /* end for i */; 508 } /* end RoundEnd */; 509 510 511 /*---------------------------------------------------------------------------- 512 * Routine: Paramaterize (xpoints, ypoints, hparams, num_points) 513 * 514 * Results: This routine calculates parameteric values for use in 515 * calculating curves. The parametric values are returned 516 * in the array h. The values are an approximation of 517 * cumulative arc lengths of the curve (uses cord length). 518 * For additional information, see paper cited below. 519 *----------------------------------------------------------------------------*/ 520 521 static Paramaterize(x, y, h, n) 522 int x[MAXPOINTS]; 523 int y[MAXPOINTS]; 524 float h[MAXPOINTS]; 525 int n; 526 { 527 register int dx; 528 register int dy; 529 register int i; 530 register int j; 531 float u[MAXPOINTS]; 532 533 534 for (i=1; i<=n; ++i) { 535 u[i] = 0; 536 for (j=1; j<i; j++) { 537 dx = x[j+1] - x[j]; 538 dy = y[j+1] - y[j]; 539 u[i] += sqrt ((double) (dx * dx + dy * dy)); 540 } 541 } 542 for (i=1; i<n; ++i) h[i] = u[i+1] - u[i]; 543 } /* end Paramaterize */ 544 545 546 /*---------------------------------------------------------------------------- 547 * Routine: PeriodicSpline (h, z, dz, d2z, d3z, npoints) 548 * 549 * Results: This routine solves for the cubic polynomial to fit a 550 * spline curve to the the points specified by the list 551 * of values. The Curve generated is periodic. The algorithms 552 * for this curve are from the "Spline Curve Techniques" paper 553 * cited below. 554 *----------------------------------------------------------------------------*/ 555 556 static PeriodicSpline(h, z, dz, d2z, d3z, npoints) 557 float h[MAXPOINTS]; /* paramaterization */ 558 int z[MAXPOINTS]; /* point list */ 559 float dz[MAXPOINTS]; /* to return the 1st derivative */ 560 float d2z[MAXPOINTS], d3z[MAXPOINTS]; /* 2nd and 3rd derivatives */ 561 int npoints; /* number of valid points */ 562 { 563 float d[MAXPOINTS]; 564 float deltaz[MAXPOINTS], a[MAXPOINTS], b[MAXPOINTS]; 565 float c[MAXPOINTS], r[MAXPOINTS], s[MAXPOINTS]; 566 int i; 567 568 /* step 1 */ 569 for (i=1; i<npoints; ++i) { 570 deltaz[i] = h[i] ? ((double) (z[i+1] - z[i])) / h[i] : 0; 571 } 572 h[0] = h[npoints-1]; 573 deltaz[0] = deltaz[npoints-1]; 574 575 /* step 2 */ 576 for (i=1; i<npoints-1; ++i) { 577 d[i] = deltaz[i+1] - deltaz[i]; 578 } 579 d[0] = deltaz[1] - deltaz[0]; 580 581 /* step 3a */ 582 a[1] = 2 * (h[0] + h[1]); 583 b[1] = d[0]; 584 c[1] = h[0]; 585 for (i=2; i<npoints-1; ++i) { 586 a[i] = 2*(h[i-1]+h[i]) - pow ((double) h[i-1],(double)2.0) / a[i-1]; 587 b[i] = d[i-1] - h[i-1] * b[i-1]/a[i-1]; 588 c[i] = -h[i-1] * c[i-1]/a[i-1]; 589 } 590 591 /* step 3b */ 592 r[npoints-1] = 1; 593 s[npoints-1] = 0; 594 for (i=npoints-2; i>0; --i) { 595 r[i] = -(h[i] * r[i+1] + c[i])/a[i]; 596 s[i] = (6 * b[i] - h[i] * s[i+1])/a[i]; 597 } 598 599 /* step 4 */ 600 d2z[npoints-1] = (6 * d[npoints-2] - h[0] * s[1] 601 - h[npoints-1] * s[npoints-2]) 602 / (h[0] * r[1] + h[npoints-1] * r[npoints-2] 603 + 2 * (h[npoints-2] + h[0])); 604 for (i=1; i<npoints-1; ++i) { 605 d2z[i] = r[i] * d2z[npoints-1] + s[i]; 606 } 607 d2z[npoints] = d2z[1]; 608 609 /* step 5 */ 610 for (i=1; i<npoints; ++i) { 611 dz[i] = deltaz[i] - h[i] * (2 * d2z[i] + d2z[i+1])/6; 612 d3z[i] = h[i] ? (d2z[i+1] - d2z[i])/h[i] : 0; 613 } 614 } /* end PeriodicSpline */ 615 616 617 /*---------------------------------------------------------------------------- 618 * Routine: NaturalEndSpline (h, z, dz, d2z, d3z, npoints) 619 * 620 * Results: This routine solves for the cubic polynomial to fit a 621 * spline curve the the points specified by the list of 622 * values. The alogrithms for this curve are from the 623 * "Spline Curve Techniques" paper cited below. 624 *----------------------------------------------------------------------------*/ 625 626 static NaturalEndSpline(h, z, dz, d2z, d3z, npoints) 627 float h[MAXPOINTS]; /* parameterization */ 628 int z[MAXPOINTS]; /* Point list */ 629 float dz[MAXPOINTS]; /* to return the 1st derivative */ 630 float d2z[MAXPOINTS], d3z[MAXPOINTS]; /* 2nd and 3rd derivatives */ 631 int npoints; /* number of valid points */ 632 { 633 float d[MAXPOINTS]; 634 float deltaz[MAXPOINTS], a[MAXPOINTS], b[MAXPOINTS]; 635 int i; 636 637 /* step 1 */ 638 for (i=1; i<npoints; ++i) { 639 deltaz[i] = h[i] ? ((double) (z[i+1] - z[i])) / h[i] : 0; 640 } 641 deltaz[0] = deltaz[npoints-1]; 642 643 /* step 2 */ 644 for (i=1; i<npoints-1; ++i) { 645 d[i] = deltaz[i+1] - deltaz[i]; 646 } 647 d[0] = deltaz[1] - deltaz[0]; 648 649 /* step 3 */ 650 a[0] = 2 * (h[2] + h[1]); 651 b[0] = d[1]; 652 for (i=1; i<npoints-2; ++i) { 653 a[i] = 2*(h[i+1]+h[i+2]) - pow((double) h[i+1],(double) 2.0)/a[i-1]; 654 b[i] = d[i+1] - h[i+1] * b[i-1]/a[i-1]; 655 } 656 657 /* step 4 */ 658 d2z[npoints] = d2z[1] = 0; 659 for (i=npoints-1; i>1; --i) { 660 d2z[i] = (6 * b[i-2] - h[i] *d2z[i+1])/a[i-2]; 661 } 662 663 /* step 5 */ 664 for (i=1; i<npoints; ++i) { 665 dz[i] = deltaz[i] - h[i] * (2 * d2z[i] + d2z[i+1])/6; 666 d3z[i] = h[i] ? (d2z[i+1] - d2z[i])/h[i] : 0; 667 } 668 } /* end NaturalEndSpline */ 669 670 671 /*---------------------------------------------------------------------------- 672 * Routine: HGCurve(xpoints, ypoints, num_points) 673 * 674 * Results: This routine generates a smooth curve through a set of points. 675 * The method used is the parametric spline curve on unit knot 676 * mesh described in "Spline Curve Techniques" by Patrick 677 * Baudelaire, Robert Flegal, and Robert Sproull -- Xerox Parc. 678 *----------------------------------------------------------------------------*/ 679 680 #define PointsPerInterval 32 681 682 HGCurve(x, y, numpoints) 683 int x[MAXPOINTS]; 684 int y[MAXPOINTS]; 685 int numpoints; 686 { 687 float h[MAXPOINTS], dx[MAXPOINTS], dy[MAXPOINTS]; 688 float d2x[MAXPOINTS], d2y[MAXPOINTS], d3x[MAXPOINTS], d3y[MAXPOINTS]; 689 float t, t2, t3; 690 register int j; 691 register int k; 692 register int nx; 693 register int ny; 694 int lx, ly; 695 696 697 lx = x[1]; 698 ly = y[1]; 699 700 /* Solve for derivatives of the curve at each point 701 * separately for x and y (parametric). 702 */ 703 Paramaterize(x, y, h, numpoints); 704 /* closed curve */ 705 if ((x[1] == x[numpoints]) && (y[1] == y[numpoints])) { 706 PeriodicSpline(h, x, dx, d2x, d3x, numpoints); 707 PeriodicSpline(h, y, dy, d2y, d3y, numpoints); 708 } else { 709 NaturalEndSpline(h, x, dx, d2x, d3x, numpoints); 710 NaturalEndSpline(h, y, dy, d2y, d3y, numpoints); 711 } 712 713 /* generate the curve using the above information and 714 * PointsPerInterval vectors between each specified knot. 715 */ 716 717 for (j=1; j<numpoints; ++j) { 718 if ((x[j] == x[j+1]) && (y[j] == y[j+1])) continue; 719 for (k=0; k<=PointsPerInterval; ++k) { 720 t = (float) k * h[j] / (float) PointsPerInterval; 721 t2 = t * t; 722 t3 = t * t * t; 723 nx = x[j] + (int) (t * dx[j] + t2 * d2x[j]/2 + t3 * d3x[j]/6); 724 ny = y[j] + (int) (t * dy[j] + t2 * d2y[j]/2 + t3 * d3y[j]/6); 725 HGtline(lx, ly, nx, ny); 726 lx = nx; 727 ly = ny; 728 } /* end for k */ 729 } /* end for j */ 730 } /* end HGCurve */ 731 732 733 /*---------------------------------------------------------------------------- 734 * Routine: HGtline(xstart, ystart, xend, yend) 735 * 736 * Results: Draws a line of proper thickness by calling "line" numerous 737 * times until the desired thickness is reached. 738 *----------------------------------------------------------------------------*/ 739 740 HGtline(x0, y0, x1, y1) 741 register int x0; 742 register int y0; 743 int x1; 744 int y1; 745 { 746 register int xs; 747 register int xe; 748 register int ys; 749 register int ye; 750 double morelen, theta, wx, wy, xx, xy; 751 int addln, j, xdir, ydir, dx, dy; 752 753 754 xdir = ydir = 1; 755 dx = x1 - x0; /* calculate direction to move to */ 756 dy = y1 - y0; /* move to draw additional lines if needed */ 757 if (dx < 0 ) { /* for extra thickness */ 758 dx = -dx; 759 xdir = -1; 760 } 761 if (dy < 0 ) { 762 dy = -dy; 763 ydir = -1; 764 } 765 766 morelen = linethickness / 2; 767 addln = (int) morelen; 768 RoundEnd (x0, y0, (int) morelen, TRUE); /* add rounded end */ 769 770 for (j=(-addln); j<=addln; ++j) { 771 if (dy == 0) { 772 xs = x0; 773 xe = x1; 774 ys = ye = y0 + j; 775 } 776 if (dx == 0) { 777 ys = y0; 778 ye = y1; 779 xs = xe = x0 + j; 780 } 781 if ((dx != 0) && (dy != 0)) { 782 theta = pi / 2.0 - atan( ((double) dx)/((double) dy) ); 783 wx = j * sin(theta); 784 wy = j * cos(theta); 785 xs = x0 + wx * xdir; 786 ys = y0 - wy * ydir; 787 xe = x1 + wx * xdir; 788 ye = y1 - wy * ydir; 789 } 790 line(xs, ys, xe, ye); 791 } /* end for */ 792 793 RoundEnd(x1, y1, (int) morelen, TRUE); /* add rounded end */ 794 } /* end HGtline */ 795