1 /* $NetBSD: hgraph.cpp,v 1.1.1.1 2016/01/13 18:41:49 christos Exp $ */ 2 3 /* Last non-groff version: hgraph.c 1.14 (Berkeley) 84/11/27 4 * 5 * This file contains the graphics routines for converting gremlin pictures 6 * to troff input. 7 */ 8 9 #include "lib.h" 10 11 #include "gprint.h" 12 13 #define MAXVECT 40 14 #define MAXPOINTS 200 15 #define LINELENGTH 1 16 #define PointsPerInterval 64 17 #define pi 3.14159265358979324 18 #define twopi (2.0 * pi) 19 #define len(a, b) groff_hypot((double)(b.x-a.x), (double)(b.y-a.y)) 20 21 22 extern int dotshifter; /* for the length of dotted curves */ 23 24 extern int style[]; /* line and character styles */ 25 extern double thick[]; 26 extern char *tfont[]; 27 extern int tsize[]; 28 extern int stipple_index[]; /* stipple font index for stipples 0 - 16 */ 29 extern char *stipple; /* stipple type (cf or ug) */ 30 31 32 extern double troffscale; /* imports from main.c */ 33 extern double linethickness; 34 extern int linmod; 35 extern int lastx; 36 extern int lasty; 37 extern int lastyline; 38 extern int ytop; 39 extern int ybottom; 40 extern int xleft; 41 extern int xright; 42 extern enum E { 43 OUTLINE, FILL, BOTH 44 } polyfill; 45 46 extern double adj1; 47 extern double adj2; 48 extern double adj3; 49 extern double adj4; 50 extern int res; 51 52 void HGSetFont(int font, int size); 53 void HGPutText(int justify, POINT pnt, register char *string); 54 void HGSetBrush(int mode); 55 void tmove2(int px, int py); 56 void doarc(POINT cp, POINT sp, int angle); 57 void tmove(POINT * ptr); 58 void cr(); 59 void drawwig(POINT * ptr, int type); 60 void HGtline(int x1, int y1); 61 void deltax(double x); 62 void deltay(double y); 63 void HGArc(register int cx, register int cy, int px, int py, int angle); 64 void picurve(register int *x, register int *y, int npts); 65 void HGCurve(int *x, int *y, int numpoints); 66 void Paramaterize(int x[], int y[], double h[], int n); 67 void PeriodicSpline(double h[], int z[], 68 double dz[], double d2z[], double d3z[], 69 int npoints); 70 void NaturalEndSpline(double h[], int z[], 71 double dz[], double d2z[], double d3z[], 72 int npoints); 73 74 75 76 /*----------------------------------------------------------------------------* 77 | Routine: HGPrintElt (element_pointer, baseline) 78 | 79 | Results: Examines a picture element and calls the appropriate 80 | routine(s) to print them according to their type. After the 81 | picture is drawn, current position is (lastx, lasty). 82 *----------------------------------------------------------------------------*/ 83 84 void 85 HGPrintElt(ELT *element, 86 int /* baseline */) 87 { 88 register POINT *p1; 89 register POINT *p2; 90 register int length; 91 register int graylevel; 92 93 if (!DBNullelt(element) && !Nullpoint((p1 = element->ptlist))) { 94 /* p1 always has first point */ 95 if (TEXT(element->type)) { 96 HGSetFont(element->brushf, element->size); 97 switch (element->size) { 98 case 1: 99 p1->y += adj1; 100 break; 101 case 2: 102 p1->y += adj2; 103 break; 104 case 3: 105 p1->y += adj3; 106 break; 107 case 4: 108 p1->y += adj4; 109 break; 110 default: 111 break; 112 } 113 HGPutText(element->type, *p1, element->textpt); 114 } else { 115 if (element->brushf) /* if there is a brush, the */ 116 HGSetBrush(element->brushf); /* graphics need it set */ 117 118 switch (element->type) { 119 120 case ARC: 121 p2 = PTNextPoint(p1); 122 tmove(p2); 123 doarc(*p1, *p2, element->size); 124 cr(); 125 break; 126 127 case CURVE: 128 length = 0; /* keep track of line length */ 129 drawwig(p1, CURVE); 130 cr(); 131 break; 132 133 case BSPLINE: 134 length = 0; /* keep track of line length */ 135 drawwig(p1, BSPLINE); 136 cr(); 137 break; 138 139 case VECTOR: 140 length = 0; /* keep track of line length so */ 141 tmove(p1); /* single lines don't get long */ 142 while (!Nullpoint((p1 = PTNextPoint(p1)))) { 143 HGtline((int) (p1->x * troffscale), 144 (int) (p1->y * troffscale)); 145 if (length++ > LINELENGTH) { 146 length = 0; 147 printf("\\\n"); 148 } 149 } /* end while */ 150 cr(); 151 break; 152 153 case POLYGON: 154 { 155 /* brushf = style of outline; size = color of fill: 156 * on first pass (polyfill=FILL), do the interior using 'P' 157 * unless size=0 158 * on second pass (polyfill=OUTLINE), do the outline using a series 159 * of vectors. It might make more sense to use \D'p ...', 160 * but there is no uniform way to specify a 'fill character' 161 * that prints as 'no fill' on all output devices (and 162 * stipple fonts). 163 * If polyfill=BOTH, just use the \D'p ...' command. 164 */ 165 double firstx = p1->x; 166 double firsty = p1->y; 167 168 length = 0; /* keep track of line length so */ 169 /* single lines don't get long */ 170 171 if (polyfill == FILL || polyfill == BOTH) { 172 /* do the interior */ 173 char command = (polyfill == BOTH && element->brushf) ? 'p' : 'P'; 174 175 /* include outline, if there is one and */ 176 /* the -p flag was set */ 177 178 /* switch based on what gremlin gives */ 179 switch (element->size) { 180 case 1: 181 graylevel = 1; 182 break; 183 case 3: 184 graylevel = 2; 185 break; 186 case 12: 187 graylevel = 3; 188 break; 189 case 14: 190 graylevel = 4; 191 break; 192 case 16: 193 graylevel = 5; 194 break; 195 case 19: 196 graylevel = 6; 197 break; 198 case 21: 199 graylevel = 7; 200 break; 201 case 23: 202 graylevel = 8; 203 break; 204 default: /* who's giving something else? */ 205 graylevel = NSTIPPLES; 206 break; 207 } 208 /* int graylevel = element->size; */ 209 210 if (graylevel < 0) 211 break; 212 if (graylevel > NSTIPPLES) 213 graylevel = NSTIPPLES; 214 printf("\\D'Fg %.3f'", 215 double(1000 - stipple_index[graylevel]) / 1000.0); 216 cr(); 217 tmove(p1); 218 printf("\\D'%c", command); 219 220 while (!Nullpoint((PTNextPoint(p1)))) { 221 p1 = PTNextPoint(p1); 222 deltax((double) p1->x); 223 deltay((double) p1->y); 224 if (length++ > LINELENGTH) { 225 length = 0; 226 printf("\\\n"); 227 } 228 } /* end while */ 229 230 /* close polygon if not done so by user */ 231 if ((firstx != p1->x) || (firsty != p1->y)) { 232 deltax((double) firstx); 233 deltay((double) firsty); 234 } 235 putchar('\''); 236 cr(); 237 break; 238 } 239 /* else polyfill == OUTLINE; only draw the outline */ 240 if (!(element->brushf)) 241 break; 242 length = 0; /* keep track of line length */ 243 tmove(p1); 244 245 while (!Nullpoint((PTNextPoint(p1)))) { 246 p1 = PTNextPoint(p1); 247 HGtline((int) (p1->x * troffscale), 248 (int) (p1->y * troffscale)); 249 if (length++ > LINELENGTH) { 250 length = 0; 251 printf("\\\n"); 252 } 253 } /* end while */ 254 255 /* close polygon if not done so by user */ 256 if ((firstx != p1->x) || (firsty != p1->y)) { 257 HGtline((int) (firstx * troffscale), 258 (int) (firsty * troffscale)); 259 } 260 cr(); 261 break; 262 } /* end case POLYGON */ 263 } /* end switch */ 264 } /* end else Text */ 265 } /* end if */ 266 } /* end PrintElt */ 267 268 269 /*----------------------------------------------------------------------------* 270 | Routine: HGPutText (justification, position_point, string) 271 | 272 | Results: Given the justification, a point to position with, and a 273 | string to put, HGPutText first sends the string into a 274 | diversion, moves to the positioning point, then outputs 275 | local vertical and horizontal motions as needed to justify 276 | the text. After all motions are done, the diversion is 277 | printed out. 278 *----------------------------------------------------------------------------*/ 279 280 void 281 HGPutText(int justify, 282 POINT pnt, 283 register char *string) 284 { 285 int savelasty = lasty; /* vertical motion for text is to be */ 286 /* ignored. Save current y here */ 287 288 printf(".nr g8 \\n(.d\n"); /* save current vertical position. */ 289 printf(".ds g9 \""); /* define string containing the text. */ 290 while (*string) { /* put out the string */ 291 if (*string == '\\' && 292 *(string + 1) == '\\') { /* one character at a */ 293 printf("\\\\\\"); /* time replacing // */ 294 string++; /* by //// to prevent */ 295 } /* interpretation at */ 296 printf("%c", *(string++)); /* printout time */ 297 } 298 printf("\n"); 299 300 tmove(&pnt); /* move to positioning point */ 301 302 switch (justify) { 303 /* local vertical motions */ 304 /* (the numbers here are used to be somewhat compatible with gprint) */ 305 case CENTLEFT: 306 case CENTCENT: 307 case CENTRIGHT: 308 printf("\\v'0.85n'"); /* down half */ 309 break; 310 311 case TOPLEFT: 312 case TOPCENT: 313 case TOPRIGHT: 314 printf("\\v'1.7n'"); /* down whole */ 315 } 316 317 switch (justify) { 318 /* local horizontal motions */ 319 case BOTCENT: 320 case CENTCENT: 321 case TOPCENT: 322 printf("\\h'-\\w'\\*(g9'u/2u'"); /* back half */ 323 break; 324 325 case BOTRIGHT: 326 case CENTRIGHT: 327 case TOPRIGHT: 328 printf("\\h'-\\w'\\*(g9'u'"); /* back whole */ 329 } 330 331 printf("\\&\\*(g9\n"); /* now print the text. */ 332 printf(".sp |\\n(g8u\n"); /* restore vertical position */ 333 lasty = savelasty; /* vertical position restored to where it */ 334 lastx = xleft; /* was before text, also horizontal is at */ 335 /* left */ 336 } /* end HGPutText */ 337 338 339 /*----------------------------------------------------------------------------* 340 | Routine: doarc (center_point, start_point, angle) 341 | 342 | Results: Produces either drawarc command or a drawcircle command 343 | depending on the angle needed to draw through. 344 *----------------------------------------------------------------------------*/ 345 346 void 347 doarc(POINT cp, 348 POINT sp, 349 int angle) 350 { 351 if (angle) /* arc with angle */ 352 HGArc((int) (cp.x * troffscale), (int) (cp.y * troffscale), 353 (int) (sp.x * troffscale), (int) (sp.y * troffscale), angle); 354 else /* a full circle (angle == 0) */ 355 HGArc((int) (cp.x * troffscale), (int) (cp.y * troffscale), 356 (int) (sp.x * troffscale), (int) (sp.y * troffscale), 0); 357 } 358 359 360 /*----------------------------------------------------------------------------* 361 | Routine: HGSetFont (font_number, Point_size) 362 | 363 | Results: ALWAYS outputs a .ft and .ps directive to troff. This is 364 | done because someone may change stuff inside a text string. 365 | Changes thickness back to default thickness. Default 366 | thickness depends on font and pointsize. 367 *----------------------------------------------------------------------------*/ 368 369 void 370 HGSetFont(int font, 371 int size) 372 { 373 printf(".ft %s\n" 374 ".ps %d\n", tfont[font - 1], tsize[size - 1]); 375 linethickness = DEFTHICK; 376 } 377 378 379 /*----------------------------------------------------------------------------* 380 | Routine: HGSetBrush (line_mode) 381 | 382 | Results: Generates the troff commands to set up the line width and 383 | style of subsequent lines. Does nothing if no change is 384 | needed. 385 | 386 | Side Efct: Sets `linmode' and `linethicknes'. 387 *----------------------------------------------------------------------------*/ 388 389 void 390 HGSetBrush(int mode) 391 { 392 register int printed = 0; 393 394 if (linmod != style[--mode]) { 395 /* Groff doesn't understand \Ds, so we take it out */ 396 /* printf ("\\D's %du'", linmod = style[mode]); */ 397 linmod = style[mode]; 398 printed = 1; 399 } 400 if (linethickness != thick[mode]) { 401 linethickness = thick[mode]; 402 printf("\\h'-%.2fp'\\D't %.2fp'", linethickness, linethickness); 403 printed = 1; 404 } 405 if (printed) 406 cr(); 407 } 408 409 410 /*----------------------------------------------------------------------------* 411 | Routine: deltax (x_destination) 412 | 413 | Results: Scales and outputs a number for delta x (with a leading 414 | space) given `lastx' and x_destination. 415 | 416 | Side Efct: Resets `lastx' to x_destination. 417 *----------------------------------------------------------------------------*/ 418 419 void 420 deltax(double x) 421 { 422 register int ix = (int) (x * troffscale); 423 424 printf(" %du", ix - lastx); 425 lastx = ix; 426 } 427 428 429 /*----------------------------------------------------------------------------* 430 | Routine: deltay (y_destination) 431 | 432 | Results: Scales and outputs a number for delta y (with a leading 433 | space) given `lastyline' and y_destination. 434 | 435 | Side Efct: Resets `lastyline' to y_destination. Since `line' vertical 436 | motions don't affect `page' ones, `lasty' isn't updated. 437 *----------------------------------------------------------------------------*/ 438 439 void 440 deltay(double y) 441 { 442 register int iy = (int) (y * troffscale); 443 444 printf(" %du", iy - lastyline); 445 lastyline = iy; 446 } 447 448 449 /*----------------------------------------------------------------------------* 450 | Routine: tmove2 (px, py) 451 | 452 | Results: Produces horizontal and vertical moves for troff given the 453 | pair of points to move to and knowing the current position. 454 | Also puts out a horizontal move to start the line. This is 455 | a variation without the .sp command. 456 *----------------------------------------------------------------------------*/ 457 458 void 459 tmove2(int px, 460 int py) 461 { 462 register int dx; 463 register int dy; 464 465 if ((dy = py - lasty)) { 466 printf("\\v'%du'", dy); 467 } 468 lastyline = lasty = py; /* lasty is always set to current */ 469 if ((dx = px - lastx)) { 470 printf("\\h'%du'", dx); 471 lastx = px; 472 } 473 } 474 475 476 /*----------------------------------------------------------------------------* 477 | Routine: tmove (point_pointer) 478 | 479 | Results: Produces horizontal and vertical moves for troff given the 480 | pointer of a point to move to and knowing the current 481 | position. Also puts out a horizontal move to start the 482 | line. 483 *----------------------------------------------------------------------------*/ 484 485 void 486 tmove(POINT * ptr) 487 { 488 register int ix = (int) (ptr->x * troffscale); 489 register int iy = (int) (ptr->y * troffscale); 490 register int dx; 491 register int dy; 492 493 if ((dy = iy - lasty)) { 494 printf(".sp %du\n", dy); 495 } 496 lastyline = lasty = iy; /* lasty is always set to current */ 497 if ((dx = ix - lastx)) { 498 printf("\\h'%du'", dx); 499 lastx = ix; 500 } 501 } 502 503 504 /*----------------------------------------------------------------------------* 505 | Routine: cr ( ) 506 | 507 | Results: Ends off an input line. `.sp -1' is also added to counteract 508 | the vertical move done at the end of text lines. 509 | 510 | Side Efct: Sets `lastx' to `xleft' for troff's return to left margin. 511 *----------------------------------------------------------------------------*/ 512 513 void 514 cr() 515 { 516 printf("\n.sp -1\n"); 517 lastx = xleft; 518 } 519 520 521 /*----------------------------------------------------------------------------* 522 | Routine: line ( ) 523 | 524 | Results: Draws a single solid line to (x,y). 525 *----------------------------------------------------------------------------*/ 526 527 void 528 line(int px, 529 int py) 530 { 531 printf("\\D'l"); 532 printf(" %du", px - lastx); 533 printf(" %du'", py - lastyline); 534 lastx = px; 535 lastyline = lasty = py; 536 } 537 538 539 /*---------------------------------------------------------------------------- 540 | Routine: drawwig (ptr, type) 541 | 542 | Results: The point sequence found in the structure pointed by ptr is 543 | placed in integer arrays for further manipulation by the 544 | existing routing. With the corresponding type parameter, 545 | either picurve or HGCurve are called. 546 *----------------------------------------------------------------------------*/ 547 548 void 549 drawwig(POINT * ptr, 550 int type) 551 { 552 register int npts; /* point list index */ 553 int x[MAXPOINTS], y[MAXPOINTS]; /* point list */ 554 555 for (npts = 1; !Nullpoint(ptr); ptr = PTNextPoint(ptr), npts++) { 556 x[npts] = (int) (ptr->x * troffscale); 557 y[npts] = (int) (ptr->y * troffscale); 558 } 559 if (--npts) { 560 if (type == CURVE) /* Use the 2 different types of curves */ 561 HGCurve(&x[0], &y[0], npts); 562 else 563 picurve(&x[0], &y[0], npts); 564 } 565 } 566 567 568 /*---------------------------------------------------------------------------- 569 | Routine: HGArc (xcenter, ycenter, xstart, ystart, angle) 570 | 571 | Results: This routine plots an arc centered about (cx, cy) counter 572 | clockwise starting from the point (px, py) through `angle' 573 | degrees. If angle is 0, a full circle is drawn. It does so 574 | by creating a draw-path around the arc whose density of 575 | points depends on the size of the arc. 576 *----------------------------------------------------------------------------*/ 577 578 void 579 HGArc(register int cx, 580 register int cy, 581 int px, 582 int py, 583 int angle) 584 { 585 double xs, ys, resolution, fullcircle; 586 int m; 587 register int mask; 588 register int extent; 589 register int nx; 590 register int ny; 591 register int length; 592 register double epsilon; 593 594 xs = px - cx; 595 ys = py - cy; 596 597 length = 0; 598 599 resolution = (1.0 + groff_hypot(xs, ys) / res) * PointsPerInterval; 600 /* mask = (1 << (int) log10(resolution + 1.0)) - 1; */ 601 (void) frexp(resolution, &m); /* A bit more elegant than log10 */ 602 for (mask = 1; mask < m; mask = mask << 1); 603 mask -= 1; 604 605 epsilon = 1.0 / resolution; 606 fullcircle = (2.0 * pi) * resolution; 607 if (angle == 0) 608 extent = (int) fullcircle; 609 else 610 extent = (int) (angle * fullcircle / 360.0); 611 612 HGtline(px, py); 613 while (--extent >= 0) { 614 xs += epsilon * ys; 615 nx = cx + (int) (xs + 0.5); 616 ys -= epsilon * xs; 617 ny = cy + (int) (ys + 0.5); 618 if (!(extent & mask)) { 619 HGtline(nx, ny); /* put out a point on circle */ 620 if (length++ > LINELENGTH) { 621 length = 0; 622 printf("\\\n"); 623 } 624 } 625 } /* end for */ 626 } /* end HGArc */ 627 628 629 /*---------------------------------------------------------------------------- 630 | Routine: picurve (xpoints, ypoints, num_of_points) 631 | 632 | Results: Draws a curve delimited by (not through) the line segments 633 | traced by (xpoints, ypoints) point list. This is the `Pic' 634 | style curve. 635 *----------------------------------------------------------------------------*/ 636 637 void 638 picurve(register int *x, 639 register int *y, 640 int npts) 641 { 642 register int nseg; /* effective resolution for each curve */ 643 register int xp; /* current point (and temporary) */ 644 register int yp; 645 int pxp, pyp; /* previous point (to make lines from) */ 646 int i; /* inner curve segment traverser */ 647 int length = 0; 648 double w; /* position factor */ 649 double t1, t2, t3; /* calculation temps */ 650 651 if (x[1] == x[npts] && y[1] == y[npts]) { 652 x[0] = x[npts - 1]; /* if the lines' ends meet, make */ 653 y[0] = y[npts - 1]; /* sure the curve meets */ 654 x[npts + 1] = x[2]; 655 y[npts + 1] = y[2]; 656 } else { /* otherwise, make the ends of the */ 657 x[0] = x[1]; /* curve touch the ending points of */ 658 y[0] = y[1]; /* the line segments */ 659 x[npts + 1] = x[npts]; 660 y[npts + 1] = y[npts]; 661 } 662 663 pxp = (x[0] + x[1]) / 2; /* make the last point pointers */ 664 pyp = (y[0] + y[1]) / 2; /* point to the start of the 1st line */ 665 tmove2(pxp, pyp); 666 667 for (; npts--; x++, y++) { /* traverse the line segments */ 668 xp = x[0] - x[1]; 669 yp = y[0] - y[1]; 670 nseg = (int) groff_hypot((double) xp, (double) yp); 671 xp = x[1] - x[2]; 672 yp = y[1] - y[2]; 673 /* `nseg' is the number of line */ 674 /* segments that will be drawn for */ 675 /* each curve segment. */ 676 nseg = (int) ((double) (nseg + (int) groff_hypot((double) xp, (double) yp)) / 677 res * PointsPerInterval); 678 679 for (i = 1; i < nseg; i++) { 680 w = (double) i / (double) nseg; 681 t1 = w * w; 682 t3 = t1 + 1.0 - (w + w); 683 t2 = 2.0 - (t3 + t1); 684 xp = (((int) (t1 * x[2] + t2 * x[1] + t3 * x[0])) + 1) / 2; 685 yp = (((int) (t1 * y[2] + t2 * y[1] + t3 * y[0])) + 1) / 2; 686 687 HGtline(xp, yp); 688 if (length++ > LINELENGTH) { 689 length = 0; 690 printf("\\\n"); 691 } 692 } 693 } 694 } 695 696 697 /*---------------------------------------------------------------------------- 698 | Routine: HGCurve(xpoints, ypoints, num_points) 699 | 700 | Results: This routine generates a smooth curve through a set of 701 | points. The method used is the parametric spline curve on 702 | unit knot mesh described in `Spline Curve Techniques' by 703 | Patrick Baudelaire, Robert Flegal, and Robert Sproull -- 704 | Xerox Parc. 705 *----------------------------------------------------------------------------*/ 706 707 void 708 HGCurve(int *x, 709 int *y, 710 int numpoints) 711 { 712 double h[MAXPOINTS], dx[MAXPOINTS], dy[MAXPOINTS]; 713 double d2x[MAXPOINTS], d2y[MAXPOINTS], d3x[MAXPOINTS], d3y[MAXPOINTS]; 714 double t, t2, t3; 715 register int j; 716 register int k; 717 register int nx; 718 register int ny; 719 int lx, ly; 720 int length = 0; 721 722 lx = x[1]; 723 ly = y[1]; 724 tmove2(lx, ly); 725 726 /* 727 * Solve for derivatives of the curve at each point separately for x and y 728 * (parametric). 729 */ 730 Paramaterize(x, y, h, numpoints); 731 732 /* closed curve */ 733 if ((x[1] == x[numpoints]) && (y[1] == y[numpoints])) { 734 PeriodicSpline(h, x, dx, d2x, d3x, numpoints); 735 PeriodicSpline(h, y, dy, d2y, d3y, numpoints); 736 } else { 737 NaturalEndSpline(h, x, dx, d2x, d3x, numpoints); 738 NaturalEndSpline(h, y, dy, d2y, d3y, numpoints); 739 } 740 741 /* 742 * generate the curve using the above information and PointsPerInterval 743 * vectors between each specified knot. 744 */ 745 746 for (j = 1; j < numpoints; ++j) { 747 if ((x[j] == x[j + 1]) && (y[j] == y[j + 1])) 748 continue; 749 for (k = 0; k <= PointsPerInterval; ++k) { 750 t = (double) k *h[j] / (double) PointsPerInterval; 751 t2 = t * t; 752 t3 = t * t * t; 753 nx = x[j] + (int) (t * dx[j] + t2 * d2x[j] / 2 + t3 * d3x[j] / 6); 754 ny = y[j] + (int) (t * dy[j] + t2 * d2y[j] / 2 + t3 * d3y[j] / 6); 755 HGtline(nx, ny); 756 if (length++ > LINELENGTH) { 757 length = 0; 758 printf("\\\n"); 759 } 760 } /* end for k */ 761 } /* end for j */ 762 } /* end HGCurve */ 763 764 765 /*---------------------------------------------------------------------------- 766 | Routine: Paramaterize (xpoints, ypoints, hparams, num_points) 767 | 768 | Results: This routine calculates parameteric values for use in 769 | calculating curves. The parametric values are returned 770 | in the array h. The values are an approximation of 771 | cumulative arc lengths of the curve (uses cord length). 772 | For additional information, see paper cited below. 773 *----------------------------------------------------------------------------*/ 774 775 void 776 Paramaterize(int x[], 777 int y[], 778 double h[], 779 int n) 780 { 781 register int dx; 782 register int dy; 783 register int i; 784 register int j; 785 double u[MAXPOINTS]; 786 787 for (i = 1; i <= n; ++i) { 788 u[i] = 0; 789 for (j = 1; j < i; j++) { 790 dx = x[j + 1] - x[j]; 791 dy = y[j + 1] - y[j]; 792 /* Here was overflowing, so I changed it. */ 793 /* u[i] += sqrt ((double) (dx * dx + dy * dy)); */ 794 u[i] += groff_hypot((double) dx, (double) dy); 795 } 796 } 797 for (i = 1; i < n; ++i) 798 h[i] = u[i + 1] - u[i]; 799 } /* end Paramaterize */ 800 801 802 /*---------------------------------------------------------------------------- 803 | Routine: PeriodicSpline (h, z, dz, d2z, d3z, npoints) 804 | 805 | Results: This routine solves for the cubic polynomial to fit a spline 806 | curve to the the points specified by the list of values. 807 | The Curve generated is periodic. The algorithms for this 808 | curve are from the `Spline Curve Techniques' paper cited 809 | above. 810 *----------------------------------------------------------------------------*/ 811 812 void 813 PeriodicSpline(double h[], /* paramaterization */ 814 int z[], /* point list */ 815 double dz[], /* to return the 1st derivative */ 816 double d2z[], /* 2nd derivative */ 817 double d3z[], /* 3rd derivative */ 818 int npoints) /* number of valid points */ 819 { 820 double d[MAXPOINTS]; 821 double deltaz[MAXPOINTS], a[MAXPOINTS], b[MAXPOINTS]; 822 double c[MAXPOINTS], r[MAXPOINTS], s[MAXPOINTS]; 823 int i; 824 825 /* step 1 */ 826 for (i = 1; i < npoints; ++i) { 827 deltaz[i] = h[i] ? ((double) (z[i + 1] - z[i])) / h[i] : 0; 828 } 829 h[0] = h[npoints - 1]; 830 deltaz[0] = deltaz[npoints - 1]; 831 832 /* step 2 */ 833 for (i = 1; i < npoints - 1; ++i) { 834 d[i] = deltaz[i + 1] - deltaz[i]; 835 } 836 d[0] = deltaz[1] - deltaz[0]; 837 838 /* step 3a */ 839 a[1] = 2 * (h[0] + h[1]); 840 b[1] = d[0]; 841 c[1] = h[0]; 842 for (i = 2; i < npoints - 1; ++i) { 843 a[i] = 2 * (h[i - 1] + h[i]) - 844 pow((double) h[i - 1], (double) 2.0) / a[i - 1]; 845 b[i] = d[i - 1] - h[i - 1] * b[i - 1] / a[i - 1]; 846 c[i] = -h[i - 1] * c[i - 1] / a[i - 1]; 847 } 848 849 /* step 3b */ 850 r[npoints - 1] = 1; 851 s[npoints - 1] = 0; 852 for (i = npoints - 2; i > 0; --i) { 853 r[i] = -(h[i] * r[i + 1] + c[i]) / a[i]; 854 s[i] = (6 * b[i] - h[i] * s[i + 1]) / a[i]; 855 } 856 857 /* step 4 */ 858 d2z[npoints - 1] = (6 * d[npoints - 2] - h[0] * s[1] 859 - h[npoints - 1] * s[npoints - 2]) 860 / (h[0] * r[1] + h[npoints - 1] * r[npoints - 2] 861 + 2 * (h[npoints - 2] + h[0])); 862 for (i = 1; i < npoints - 1; ++i) { 863 d2z[i] = r[i] * d2z[npoints - 1] + s[i]; 864 } 865 d2z[npoints] = d2z[1]; 866 867 /* step 5 */ 868 for (i = 1; i < npoints; ++i) { 869 dz[i] = deltaz[i] - h[i] * (2 * d2z[i] + d2z[i + 1]) / 6; 870 d3z[i] = h[i] ? (d2z[i + 1] - d2z[i]) / h[i] : 0; 871 } 872 } /* end PeriodicSpline */ 873 874 875 /*---------------------------------------------------------------------------- 876 | Routine: NaturalEndSpline (h, z, dz, d2z, d3z, npoints) 877 | 878 | Results: This routine solves for the cubic polynomial to fit a spline 879 | curve the the points specified by the list of values. The 880 | alogrithms for this curve are from the `Spline Curve 881 | Techniques' paper cited above. 882 *----------------------------------------------------------------------------*/ 883 884 void 885 NaturalEndSpline(double h[], /* parameterization */ 886 int z[], /* Point list */ 887 double dz[], /* to return the 1st derivative */ 888 double d2z[], /* 2nd derivative */ 889 double d3z[], /* 3rd derivative */ 890 int npoints) /* number of valid points */ 891 { 892 double d[MAXPOINTS]; 893 double deltaz[MAXPOINTS], a[MAXPOINTS], b[MAXPOINTS]; 894 int i; 895 896 /* step 1 */ 897 for (i = 1; i < npoints; ++i) { 898 deltaz[i] = h[i] ? ((double) (z[i + 1] - z[i])) / h[i] : 0; 899 } 900 deltaz[0] = deltaz[npoints - 1]; 901 902 /* step 2 */ 903 for (i = 1; i < npoints - 1; ++i) { 904 d[i] = deltaz[i + 1] - deltaz[i]; 905 } 906 d[0] = deltaz[1] - deltaz[0]; 907 908 /* step 3 */ 909 a[0] = 2 * (h[2] + h[1]); 910 b[0] = d[1]; 911 for (i = 1; i < npoints - 2; ++i) { 912 a[i] = 2 * (h[i + 1] + h[i + 2]) - 913 pow((double) h[i + 1], (double) 2.0) / a[i - 1]; 914 b[i] = d[i + 1] - h[i + 1] * b[i - 1] / a[i - 1]; 915 } 916 917 /* step 4 */ 918 d2z[npoints] = d2z[1] = 0; 919 for (i = npoints - 1; i > 1; --i) { 920 d2z[i] = (6 * b[i - 2] - h[i] * d2z[i + 1]) / a[i - 2]; 921 } 922 923 /* step 5 */ 924 for (i = 1; i < npoints; ++i) { 925 dz[i] = deltaz[i] - h[i] * (2 * d2z[i] + d2z[i + 1]) / 6; 926 d3z[i] = h[i] ? (d2z[i + 1] - d2z[i]) / h[i] : 0; 927 } 928 } /* end NaturalEndSpline */ 929 930 931 /*----------------------------------------------------------------------------* 932 | Routine: change (x_position, y_position, visible_flag) 933 | 934 | Results: As HGtline passes from the invisible to visible (or vice 935 | versa) portion of a line, change is called to either draw 936 | the line, or initialize the beginning of the next one. 937 | Change calls line to draw segments if visible_flag is set 938 | (which means we're leaving a visible area). 939 *----------------------------------------------------------------------------*/ 940 941 void 942 change(register int x, 943 register int y, 944 register int vis) 945 { 946 static int length = 0; 947 948 if (vis) { /* leaving a visible area, draw it. */ 949 line(x, y); 950 if (length++ > LINELENGTH) { 951 length = 0; 952 printf("\\\n"); 953 } 954 } else { /* otherwise, we're entering one, remember */ 955 /* beginning */ 956 tmove2(x, y); 957 } 958 } 959 960 961 /*---------------------------------------------------------------------------- 962 | Routine: HGtline (xstart, ystart, xend, yend) 963 | 964 | Results: Draws a line from current position to (x1,y1) using line(x1, 965 | y1) to place individual segments of dotted or dashed lines. 966 *----------------------------------------------------------------------------*/ 967 968 void 969 HGtline(int x_1, 970 int y_1) 971 { 972 register int x_0 = lastx; 973 register int y_0 = lasty; 974 register int dx; 975 register int dy; 976 register int oldcoord; 977 register int res1; 978 register int visible; 979 register int res2; 980 register int xinc; 981 register int yinc; 982 register int dotcounter; 983 984 if (linmod == SOLID) { 985 line(x_1, y_1); 986 return; 987 } 988 989 /* for handling different resolutions */ 990 dotcounter = linmod << dotshifter; 991 992 xinc = 1; 993 yinc = 1; 994 if ((dx = x_1 - x_0) < 0) { 995 xinc = -xinc; 996 dx = -dx; 997 } 998 if ((dy = y_1 - y_0) < 0) { 999 yinc = -yinc; 1000 dy = -dy; 1001 } 1002 res1 = 0; 1003 res2 = 0; 1004 visible = 0; 1005 if (dx >= dy) { 1006 oldcoord = y_0; 1007 while (x_0 != x_1) { 1008 if ((x_0 & dotcounter) && !visible) { 1009 change(x_0, y_0, 0); 1010 visible = 1; 1011 } else if (visible && !(x_0 & dotcounter)) { 1012 change(x_0 - xinc, oldcoord, 1); 1013 visible = 0; 1014 } 1015 if (res1 > res2) { 1016 oldcoord = y_0; 1017 res2 += dx - res1; 1018 res1 = 0; 1019 y_0 += yinc; 1020 } 1021 res1 += dy; 1022 x_0 += xinc; 1023 } 1024 } else { 1025 oldcoord = x_0; 1026 while (y_0 != y_1) { 1027 if ((y_0 & dotcounter) && !visible) { 1028 change(x_0, y_0, 0); 1029 visible = 1; 1030 } else if (visible && !(y_0 & dotcounter)) { 1031 change(oldcoord, y_0 - yinc, 1); 1032 visible = 0; 1033 } 1034 if (res1 > res2) { 1035 oldcoord = x_0; 1036 res2 += dy - res1; 1037 res1 = 0; 1038 x_0 += xinc; 1039 } 1040 res1 += dx; 1041 y_0 += yinc; 1042 } 1043 } 1044 if (visible) 1045 change(x_1, y_1, 1); 1046 else 1047 change(x_1, y_1, 0); 1048 } 1049 1050 /* EOF */ 1051