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