Sync usage with man page.
[netbsd-mini2440.git] / gnu / dist / groff / src / preproc / grn / hgraph.cpp
blobd81ac937a72c806655448a02cfd17cb61cc28b1d
1 /* $NetBSD$ */
3 /* Last non-groff version: hgraph.c 1.14 (Berkeley) 84/11/27
5 * This file contains the graphics routines for converting gremlin pictures
6 * to troff input.
7 */
9 #include "lib.h"
11 #include "gprint.h"
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))
22 extern int dotshifter; /* for the length of dotted curves */
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) */
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;
46 extern double adj1;
47 extern double adj2;
48 extern double adj3;
49 extern double adj4;
50 extern int res;
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);
76 /*----------------------------------------------------------------------------*
77 | Routine: HGPrintElt (element_pointer, baseline)
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 *----------------------------------------------------------------------------*/
84 void
85 HGPrintElt(ELT *element,
86 int /* baseline */)
88 register POINT *p1;
89 register POINT *p2;
90 register int length;
91 register int graylevel;
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;
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 */
118 switch (element->type) {
120 case ARC:
121 p2 = PTNextPoint(p1);
122 tmove(p2);
123 doarc(*p1, *p2, element->size);
124 cr();
125 break;
127 case CURVE:
128 length = 0; /* keep track of line length */
129 drawwig(p1, CURVE);
130 cr();
131 break;
133 case BSPLINE:
134 length = 0; /* keep track of line length */
135 drawwig(p1, BSPLINE);
136 cr();
137 break;
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");
149 } /* end while */
150 cr();
151 break;
153 case POLYGON:
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.
165 double firstx = p1->x;
166 double firsty = p1->y;
168 length = 0; /* keep track of line length so */
169 /* single lines don't get long */
171 if (polyfill == FILL || polyfill == BOTH) {
172 /* do the interior */
173 char command = (polyfill == BOTH && element->brushf) ? 'p' : 'P';
175 /* include outline, if there is one and */
176 /* the -p flag was set */
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;
208 /* int graylevel = element->size; */
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);
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");
228 } /* end while */
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);
235 putchar('\'');
236 cr();
237 break;
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);
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");
253 } /* end while */
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));
260 cr();
261 break;
262 } /* end case POLYGON */
263 } /* end switch */
264 } /* end else Text */
265 } /* end if */
266 } /* end PrintElt */
269 /*----------------------------------------------------------------------------*
270 | Routine: HGPutText (justification, position_point, string)
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 *----------------------------------------------------------------------------*/
280 void
281 HGPutText(int justify,
282 POINT pnt,
283 register char *string)
285 int savelasty = lasty; /* vertical motion for text is to be */
286 /* ignored. Save current y here */
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 */
298 printf("\n");
300 tmove(&pnt); /* move to positioning point */
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;
311 case TOPLEFT:
312 case TOPCENT:
313 case TOPRIGHT:
314 printf("\\v'1.7n'"); /* down whole */
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;
325 case BOTRIGHT:
326 case CENTRIGHT:
327 case TOPRIGHT:
328 printf("\\h'-\\w'\\*(g9'u'"); /* back whole */
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 */
339 /*----------------------------------------------------------------------------*
340 | Routine: doarc (center_point, start_point, angle)
342 | Results: Produces either drawarc command or a drawcircle command
343 | depending on the angle needed to draw through.
344 *----------------------------------------------------------------------------*/
346 void
347 doarc(POINT cp,
348 POINT sp,
349 int angle)
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);
360 /*----------------------------------------------------------------------------*
361 | Routine: HGSetFont (font_number, Point_size)
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 *----------------------------------------------------------------------------*/
369 void
370 HGSetFont(int font,
371 int size)
373 printf(".ft %s\n"
374 ".ps %d\n", tfont[font - 1], tsize[size - 1]);
375 linethickness = DEFTHICK;
379 /*----------------------------------------------------------------------------*
380 | Routine: HGSetBrush (line_mode)
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.
386 | Side Efct: Sets `linmode' and `linethicknes'.
387 *----------------------------------------------------------------------------*/
389 void
390 HGSetBrush(int mode)
392 register int printed = 0;
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;
400 if (linethickness != thick[mode]) {
401 linethickness = thick[mode];
402 printf("\\h'-%.2fp'\\D't %.2fp'", linethickness, linethickness);
403 printed = 1;
405 if (printed)
406 cr();
410 /*----------------------------------------------------------------------------*
411 | Routine: deltax (x_destination)
413 | Results: Scales and outputs a number for delta x (with a leading
414 | space) given `lastx' and x_destination.
416 | Side Efct: Resets `lastx' to x_destination.
417 *----------------------------------------------------------------------------*/
419 void
420 deltax(double x)
422 register int ix = (int) (x * troffscale);
424 printf(" %du", ix - lastx);
425 lastx = ix;
429 /*----------------------------------------------------------------------------*
430 | Routine: deltay (y_destination)
432 | Results: Scales and outputs a number for delta y (with a leading
433 | space) given `lastyline' and y_destination.
435 | Side Efct: Resets `lastyline' to y_destination. Since `line' vertical
436 | motions don't affect `page' ones, `lasty' isn't updated.
437 *----------------------------------------------------------------------------*/
439 void
440 deltay(double y)
442 register int iy = (int) (y * troffscale);
444 printf(" %du", iy - lastyline);
445 lastyline = iy;
449 /*----------------------------------------------------------------------------*
450 | Routine: tmove2 (px, py)
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 *----------------------------------------------------------------------------*/
458 void
459 tmove2(int px,
460 int py)
462 register int dx;
463 register int dy;
465 if ((dy = py - lasty)) {
466 printf("\\v'%du'", dy);
468 lastyline = lasty = py; /* lasty is always set to current */
469 if ((dx = px - lastx)) {
470 printf("\\h'%du'", dx);
471 lastx = px;
476 /*----------------------------------------------------------------------------*
477 | Routine: tmove (point_pointer)
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 *----------------------------------------------------------------------------*/
485 void
486 tmove(POINT * ptr)
488 register int ix = (int) (ptr->x * troffscale);
489 register int iy = (int) (ptr->y * troffscale);
490 register int dx;
491 register int dy;
493 if ((dy = iy - lasty)) {
494 printf(".sp %du\n", dy);
496 lastyline = lasty = iy; /* lasty is always set to current */
497 if ((dx = ix - lastx)) {
498 printf("\\h'%du'", dx);
499 lastx = ix;
504 /*----------------------------------------------------------------------------*
505 | Routine: cr ( )
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.
510 | Side Efct: Sets `lastx' to `xleft' for troff's return to left margin.
511 *----------------------------------------------------------------------------*/
513 void
514 cr()
516 printf("\n.sp -1\n");
517 lastx = xleft;
521 /*----------------------------------------------------------------------------*
522 | Routine: line ( )
524 | Results: Draws a single solid line to (x,y).
525 *----------------------------------------------------------------------------*/
527 void
528 line(int px,
529 int py)
531 printf("\\D'l");
532 printf(" %du", px - lastx);
533 printf(" %du'", py - lastyline);
534 lastx = px;
535 lastyline = lasty = py;
539 /*----------------------------------------------------------------------------
540 | Routine: drawwig (ptr, type)
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 *----------------------------------------------------------------------------*/
548 void
549 drawwig(POINT * ptr,
550 int type)
552 register int npts; /* point list index */
553 int x[MAXPOINTS], y[MAXPOINTS]; /* point list */
555 for (npts = 1; !Nullpoint(ptr); ptr = PTNextPoint(ptr), npts++) {
556 x[npts] = (int) (ptr->x * troffscale);
557 y[npts] = (int) (ptr->y * troffscale);
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);
568 /*----------------------------------------------------------------------------
569 | Routine: HGArc (xcenter, ycenter, xstart, ystart, angle)
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 *----------------------------------------------------------------------------*/
578 void
579 HGArc(register int cx,
580 register int cy,
581 int px,
582 int py,
583 int angle)
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;
594 xs = px - cx;
595 ys = py - cy;
597 length = 0;
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;
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);
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");
625 } /* end for */
626 } /* end HGArc */
629 /*----------------------------------------------------------------------------
630 | Routine: picurve (xpoints, ypoints, num_of_points)
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 *----------------------------------------------------------------------------*/
637 void
638 picurve(register int *x,
639 register int *y,
640 int npts)
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 */
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];
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);
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);
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;
687 HGtline(xp, yp);
688 if (length++ > LINELENGTH) {
689 length = 0;
690 printf("\\\n");
697 /*----------------------------------------------------------------------------
698 | Routine: HGCurve(xpoints, ypoints, num_points)
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 *----------------------------------------------------------------------------*/
707 void
708 HGCurve(int *x,
709 int *y,
710 int numpoints)
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;
722 lx = x[1];
723 ly = y[1];
724 tmove2(lx, ly);
727 * Solve for derivatives of the curve at each point separately for x and y
728 * (parametric).
730 Paramaterize(x, y, h, numpoints);
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);
742 * generate the curve using the above information and PointsPerInterval
743 * vectors between each specified knot.
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");
760 } /* end for k */
761 } /* end for j */
762 } /* end HGCurve */
765 /*----------------------------------------------------------------------------
766 | Routine: Paramaterize (xpoints, ypoints, hparams, num_points)
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 *----------------------------------------------------------------------------*/
775 void
776 Paramaterize(int x[],
777 int y[],
778 double h[],
779 int n)
781 register int dx;
782 register int dy;
783 register int i;
784 register int j;
785 double u[MAXPOINTS];
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);
797 for (i = 1; i < n; ++i)
798 h[i] = u[i + 1] - u[i];
799 } /* end Paramaterize */
802 /*----------------------------------------------------------------------------
803 | Routine: PeriodicSpline (h, z, dz, d2z, d3z, npoints)
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 *----------------------------------------------------------------------------*/
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 */
820 double d[MAXPOINTS];
821 double deltaz[MAXPOINTS], a[MAXPOINTS], b[MAXPOINTS];
822 double c[MAXPOINTS], r[MAXPOINTS], s[MAXPOINTS];
823 int i;
825 /* step 1 */
826 for (i = 1; i < npoints; ++i) {
827 deltaz[i] = h[i] ? ((double) (z[i + 1] - z[i])) / h[i] : 0;
829 h[0] = h[npoints - 1];
830 deltaz[0] = deltaz[npoints - 1];
832 /* step 2 */
833 for (i = 1; i < npoints - 1; ++i) {
834 d[i] = deltaz[i + 1] - deltaz[i];
836 d[0] = deltaz[1] - deltaz[0];
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];
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];
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];
865 d2z[npoints] = d2z[1];
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;
872 } /* end PeriodicSpline */
875 /*----------------------------------------------------------------------------
876 | Routine: NaturalEndSpline (h, z, dz, d2z, d3z, npoints)
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 *----------------------------------------------------------------------------*/
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 */
892 double d[MAXPOINTS];
893 double deltaz[MAXPOINTS], a[MAXPOINTS], b[MAXPOINTS];
894 int i;
896 /* step 1 */
897 for (i = 1; i < npoints; ++i) {
898 deltaz[i] = h[i] ? ((double) (z[i + 1] - z[i])) / h[i] : 0;
900 deltaz[0] = deltaz[npoints - 1];
902 /* step 2 */
903 for (i = 1; i < npoints - 1; ++i) {
904 d[i] = deltaz[i + 1] - deltaz[i];
906 d[0] = deltaz[1] - deltaz[0];
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];
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];
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;
928 } /* end NaturalEndSpline */
931 /*----------------------------------------------------------------------------*
932 | Routine: change (x_position, y_position, visible_flag)
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 *----------------------------------------------------------------------------*/
941 void
942 change(register int x,
943 register int y,
944 register int vis)
946 static int length = 0;
948 if (vis) { /* leaving a visible area, draw it. */
949 line(x, y);
950 if (length++ > LINELENGTH) {
951 length = 0;
952 printf("\\\n");
954 } else { /* otherwise, we're entering one, remember */
955 /* beginning */
956 tmove2(x, y);
961 /*----------------------------------------------------------------------------
962 | Routine: HGtline (xstart, ystart, xend, yend)
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 *----------------------------------------------------------------------------*/
968 void
969 HGtline(int x_1,
970 int y_1)
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;
984 if (linmod == SOLID) {
985 line(x_1, y_1);
986 return;
989 /* for handling different resolutions */
990 dotcounter = linmod << dotshifter;
992 xinc = 1;
993 yinc = 1;
994 if ((dx = x_1 - x_0) < 0) {
995 xinc = -xinc;
996 dx = -dx;
998 if ((dy = y_1 - y_0) < 0) {
999 yinc = -yinc;
1000 dy = -dy;
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;
1015 if (res1 > res2) {
1016 oldcoord = y_0;
1017 res2 += dx - res1;
1018 res1 = 0;
1019 y_0 += yinc;
1021 res1 += dy;
1022 x_0 += xinc;
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;
1034 if (res1 > res2) {
1035 oldcoord = x_0;
1036 res2 += dy - res1;
1037 res1 = 0;
1038 x_0 += xinc;
1040 res1 += dx;
1041 y_0 += yinc;
1044 if (visible)
1045 change(x_1, y_1, 1);
1046 else
1047 change(x_1, y_1, 0);
1050 /* EOF */