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