1 /*****************************************************************************/
2 /* */
3 /* ,d88^^o 888 o o */
4 /* 8888 888o^88, o88^^o Y88b o / d8b d8b o88^^8o */
5 /* "Y88b 888 888 d888 b Y88b d8b / d888bdY88b d888 88b */
6 /* "Y88b, 888 888 8888 8 Y888/Y88b/ / Y88Y Y888b 8888oo888 */
7 /* o 8888 888 888 q888 p Y8/ Y8/ / YY Y888b q888 */
8 /* "oo88P" 888 888 "88oo" Y Y / Y888b "88oooo" */
9 /* */
10 /* A Display Program for Meshes and More. */
11 /* (showme.c) */
12 /* */
13 /* Version 1.3 */
14 /* July 20, 1996 */
15 /* */
16 /* Copyright 1996 */
17 /* Jonathan Richard Shewchuk */
18 /* School of Computer Science */
19 /* Carnegie Mellon University */
20 /* 5000 Forbes Avenue */
21 /* Pittsburgh, Pennsylvania 15213-3891 */
22 /* jrs@cs.cmu.edu */
23 /* */
24 /* This program may be freely redistributed under the condition that the */
25 /* copyright notices (including this entire header and the copyright */
26 /* notice printed when the `-h' switch is selected) are not removed, and */
27 /* no compensation is received. Private, research, and institutional */
28 /* use is free. You may distribute modified versions of this code UNDER */
29 /* THE CONDITION THAT THIS CODE AND ANY MODIFICATIONS MADE TO IT IN THE */
30 /* SAME FILE REMAIN UNDER COPYRIGHT OF THE ORIGINAL AUTHOR, BOTH SOURCE */
31 /* AND OBJECT CODE ARE MADE FREELY AVAILABLE WITHOUT CHARGE, AND CLEAR */
32 /* NOTICE IS GIVEN OF THE MODIFICATIONS. Distribution of this code as */
33 /* part of a commercial system is permissible ONLY BY DIRECT ARRANGEMENT */
34 /* WITH THE AUTHOR. (If you are not directly supplying this code to a */
35 /* customer, and you are instead telling them how they can obtain it for */
36 /* free, then you are not required to make any arrangement with me.) */
37 /* */
38 /* Hypertext instructions for Triangle are available on the Web at */
39 /* */
40 /* http://www.cs.cmu.edu/~quake/showme.html */
41 /* */
42 /* Show Me was created as part of the Archimedes project in the School of */
43 /* Computer Science at Carnegie Mellon University. Archimedes is a */
44 /* system for compiling parallel finite element solvers. For further */
45 /* information, see Anja Feldmann, Omar Ghattas, John R. Gilbert, Gary L. */
46 /* Miller, David R. O'Hallaron, Eric J. Schwabe, Jonathan R. Shewchuk, */
47 /* and Shang-Hua Teng. "Automated Parallel Solution of Unstructured PDE */
48 /* Problems." To appear in Communications of the ACM, we hope. */
49 /* */
50 /* If you make any improvements to this code, please please please let me */
51 /* know, so that I may obtain the improvements. Even if you don't change */
52 /* the code, I'd still love to hear what it's being used for. */
53 /* */
54 /* Disclaimer: Neither I nor Carnegie Mellon warrant this code in any way */
55 /* whatsoever. Use at your own risk. */
56 /* */
57 /*****************************************************************************/
58
59 /* For single precision (which will save some memory and reduce paging), */
60 /* write "#define SINGLE" below. */
61 /* */
62 /* For double precision (which will allow you to display triangulations of */
63 /* a finer resolution), leave SINGLE undefined. */
64
65 /* #define SINGLE */
66
67 #ifdef SINGLE
68 #define REAL float
69 #else
70 #define REAL double
71 #endif
72
73 /* Maximum number of characters in a file name (including the null). */
74
75 #define FILENAMESIZE 1024
76
77 /* Maximum number of characters in a line read from a file (including the */
78 /* null). */
79
80 #define INPUTLINESIZE 512
81
82 #define STARTWIDTH 414
83 #define STARTHEIGHT 414
84 #define MINWIDTH 50
85 #define MINHEIGHT 50
86 #define BUTTONHEIGHT 21
87 #define BUTTONROWS 3
88 #define PANELHEIGHT (BUTTONHEIGHT * BUTTONROWS)
89 #define MAXCOLORS 64
90
91 #define IMAGE_TYPES 7
92 #define NOTHING -1
93 #define NODE 0
94 #define POLY 1
95 #define ELE 2
96 #define EDGE 3
97 #define PART 4
98 #define ADJ 5
99 #define VORO 6
100
101 #define STARTEXPLOSION 0.5
102
103 #include <stdio.h>
104 #include <stdlib.h>
105 #include <string.h>
106 #include <X11/Xlib.h>
107 #include <X11/Xutil.h>
108 #include <X11/Xatom.h>
109
110 /* The following obscenity seems to be necessary to ensure that this program */
111 /* will port to Dec Alphas running OSF/1, because their stdio.h file commits */
112 /* the unpardonable sin of including stdlib.h. Hence, malloc(), free(), and */
113 /* exit() may or may not already be defined at this point. I declare these */
114 /* functions explicitly because some non-ANSI C compilers lack stdlib.h. */
115
116 #if !defined(_STDLIB_H_) && !defined(_STDLIB_H) && defined(__need_malloc_and_calloc)
117 extern char *malloc();
118 extern void free();
119 extern void exit();
120 extern double strtod();
121 extern long strtol();
122 #endif
123
124 /* A necessary forward declaration. */
125
126 int load_image();
127
128 Display *display;
129 int screen;
130 Window rootwindow;
131 Window mainwindow;
132 Window quitwin;
133 Window leftwin;
134 Window rightwin;
135 Window upwin;
136 Window downwin;
137 Window resetwin;
138 Window pswin;
139 Window epswin;
140 Window expwin;
141 Window exppluswin;
142 Window expminuswin;
143 Window widthpluswin;
144 Window widthminuswin;
145 Window versionpluswin;
146 Window versionminuswin;
147 Window fillwin;
148 Window nodewin[2];
149 Window polywin[2];
150 Window elewin[2];
151 Window edgewin[2];
152 Window partwin[2];
153 Window adjwin[2];
154 Window voronoiwin[2];
155
156 int windowdepth;
157 XEvent event;
158 Colormap rootmap;
159 XFontStruct *font;
160 int width, height;
161 int black, white;
162 int showme_foreground;
163 GC fontgc;
164 GC blackfontgc;
165 GC linegc;
166 GC trianglegc;
167 int colors[MAXCOLORS];
168 XColor rgb[MAXCOLORS];
169 int color;
170
171 int start_image, current_image;
172 int start_inc, current_inc;
173 int loweriteration;
174 int line_width;
175 int loaded[2][IMAGE_TYPES];
176 REAL xlo[2][IMAGE_TYPES], ylo[2][IMAGE_TYPES];
177 REAL xhi[2][IMAGE_TYPES], yhi[2][IMAGE_TYPES];
178 REAL xscale, yscale;
179 REAL xoffset, yoffset;
180 int zoom;
181
182 int nodes[2], node_dim[2];
183 REAL *nodeptr[2];
184 int polynodes[2], poly_dim[2], polyedges[2], polyholes[2];
185 REAL *polynodeptr[2], *polyholeptr[2];
186 int *polyedgeptr[2];
187 int elems[2], ele_corners[2];
188 int *eleptr[2];
189 int edges[2];
190 int *edgeptr[2];
191 REAL *normptr[2];
192 int subdomains[2];
193 int *partpart[2];
194 REAL *partcenter[2], *partshift[2];
195 int adjsubdomains[2];
196 int *adjptr[2];
197 int vnodes[2], vnode_dim[2];
198 REAL *vnodeptr[2];
199 int vedges[2];
200 int *vedgeptr[2];
201 REAL *vnormptr[2];
202 int firstnumber[2];
203
204 int quiet, fillelem, bw_ps, explode;
205 REAL explosion;
206
207 char filename[FILENAMESIZE];
208 char nodefilename[2][FILENAMESIZE];
209 char polyfilename[2][FILENAMESIZE];
210 char elefilename[2][FILENAMESIZE];
211 char edgefilename[2][FILENAMESIZE];
212 char partfilename[2][FILENAMESIZE];
213 char adjfilename[2][FILENAMESIZE];
214 char vnodefilename[2][FILENAMESIZE];
215 char vedgefilename[2][FILENAMESIZE];
216
217 const
218 char *colorname[] = {"aquamarine", "red", "green yellow", "magenta",
219 "yellow", "green", "orange", "blue",
220 "white", "sandy brown", "cyan", "moccasin",
221 "cadet blue", "coral", "cornflower blue", "sky blue",
222 "firebrick", "forest green", "gold", "goldenrod",
223 "gray", "hot pink", "chartreuse", "pale violet red",
224 "indian red", "khaki", "lavender", "light blue",
225 "light gray", "light steel blue", "lime green", "azure",
226 "maroon", "medium aquamarine", "dodger blue", "honeydew",
227 "medium orchid", "medium sea green", "moccasin",
228 "medium slate blue", "medium spring green",
229 "medium turquoise", "medium violet red",
230 "orange red", "chocolate", "light goldenrod",
231 "orchid", "pale green", "pink", "plum",
232 "purple", "salmon", "sea green",
233 "sienna", "slate blue", "spring green",
234 "steel blue", "tan", "thistle", "turquoise",
235 "violet", "violet red", "wheat",
236 "yellow green"};
237
syntax()238 void syntax()
239 {
240 printf("showme [-bfw_Qh] input_file\n");
241 printf(" -b Black and white PostScript (default is color).\n");
242 printf(" -f Fill triangles of partitioned mesh with color.\n");
243 printf(" -w Set line width to some specified number.\n");
244 printf(" -Q Quiet: No terminal output except errors.\n");
245 printf(" -h Help: Detailed instructions for Show Me.\n");
246 exit(0);
247 }
248
info()249 void info()
250 {
251 printf("Show Me\n");
252 printf("A Display Program for Meshes and More.\n");
253 printf("Version 1.3\n\n");
254 printf(
255 "Copyright 1996 Jonathan Richard Shewchuk (bugs/comments to jrs@cs.cmu.edu)\n"
256 );
257 printf("School of Computer Science / Carnegie Mellon University\n");
258 printf("5000 Forbes Avenue / Pittsburgh, Pennsylvania 15213-3891\n");
259 printf(
260 "Created as part of the Archimedes project (tools for parallel FEM).\n");
261 printf(
262 "Supported in part by NSF Grant CMS-9318163 and an NSERC 1967 Scholarship.\n");
263 printf("There is no warranty whatsoever. Use at your own risk.\n");
264 #ifdef SINGLE
265 printf("This executable is compiled for single precision arithmetic.\n\n\n");
266 #else
267 printf("This executable is compiled for double precision arithmetic.\n\n\n");
268 #endif
269 printf(
270 "Show Me graphically displays the contents of geometric files, especially\n");
271 printf(
272 "those generated by Triangle, my two-dimensional quality mesh generator and\n"
273 );
274 printf(
275 "Delaunay triangulator. Show Me can also write images in PostScript form.\n");
276 printf(
277 "Show Me is also useful for checking the consistency of the files you create\n"
278 );
279 printf(
280 "as input to Triangle; Show Me does these checks more thoroughly than\n");
281 printf("Triangle does. The command syntax is:\n\n");
282 printf("showme [-bfw_Qh] input_file\n\n");
283 printf(
284 "The underscore indicates that a number should follow the -w switch.\n");
285 printf(
286 "input_file may be one of several types of file. It must have extension\n");
287 printf(
288 ".node, .poly, .ele, .edge, .part, or .adj. If no extension is provided,\n");
289 printf(
290 "Show Me will assume the extension .ele. A .node file represents a set of\n");
291 printf(
292 "points; a .poly file represents a Planar Straight Line Graph; an .ele file\n"
293 );
294 printf(
295 "(coupled with a .node file) represents the elements of a mesh or the\n");
296 printf(
297 "triangles of a triangulation; an .edge file (coupled with a .node file)\n");
298 printf(
299 "represents a set of edges; a .part file specifies a partition of a mesh;\n");
300 printf(
301 "and a .adj file represents the adjacency graph defined by a partition.\n");
302 printf("\n");
303 printf("Command Line Switches:\n");
304 printf("\n");
305 printf(
306 " -b Makes all PostScript output black and white. If this switch is not\n"
307 );
308 printf(
309 " selected, color PostScript is used for partitioned meshes and\n");
310 printf(" adjacency graphs (.part and .adj files).\n");
311 printf(
312 " -f On color displays and in color PostScript, displays partitioned\n");
313 printf(
314 " meshes by filling triangles with color, rather than by coloring the\n"
315 );
316 printf(
317 " edges. This switch will result in a clearer picture if all\n");
318 printf(
319 " triangles are reasonably large, and a less clear picture if small\n");
320 printf(
321 " triangles are present. (There is also a button to toggle this\n");
322 printf(" behavior.)\n");
323 printf(
324 " -w Followed by an integer, specifies the line width used in all\n");
325 printf(
326 " images. (There are also buttons to change the line width.)\n");
327 printf(
328 " -Q Quiet: Suppresses all explanation of what Show Me is doing, unless\n"
329 );
330 printf(" an error occurs.\n");
331 printf(" -h Help: Displays these instructions.\n");
332 printf("\n");
333 printf("Controls:\n");
334 printf("\n");
335 printf(
336 " To zoom in on an image, point at the location where you want a closer\n");
337 printf(
338 " look, and click the left mouse button. To zoom out, click the right\n");
339 printf(
340 " mouse button. In either case, the point you click on will be centered in\n"
341 );
342 printf(
343 " the window. If you want to know the coordinates of a point, click the\n");
344 printf(
345 " middle mouse button; the coordinates will be printed on the terminal you\n"
346 );
347 printf(" invoked Show Me from.\n\n");
348 printf(
349 " If you resize the window, the image will grow or shrink to match.\n");
350 printf("\n");
351 printf(
352 " There is a panel of control buttons at the bottom of the Show Me window:\n"
353 );
354 printf("\n");
355 printf(" Quit: Shuts down Show Me.\n");
356 printf(" <, >, ^, v: Moves the image in the indicated direction.\n");
357 printf(
358 " Reset: Unzooms and centers the image in the window. When you switch from\n"
359 );
360 printf(
361 " one image to another, the viewing region does not change, so you may\n");
362 printf(
363 " need to reset the new image to make it fully visible. This often is\n");
364 printf(
365 " the case when switching between Delaunay triangulations and their\n");
366 printf(
367 " corresponding Voronoi diagrams, as Voronoi vertices can be far from the\n"
368 );
369 printf(" initial point set.\n");
370 printf(
371 " Width+, -: Increases or decreases the width of all lines and points.\n");
372 printf(
373 " Exp, +, -: These buttons appear only when you are viewing a partitioned\n"
374 );
375 printf(
376 " mesh (.part file). `Exp' toggles between an exploded and non-exploded\n"
377 );
378 printf(
379 " image of the mesh. The non-exploded image will not show the partition\n"
380 );
381 printf(
382 " on a black and white monitor. `+' and `-' allow you to adjust the\n");
383 printf(
384 " spacing between pieces of the mesh to better distinguish them.\n");
385 printf(
386 " Fill: This button appears only when you are viewing a partitioned mesh\n");
387 printf(
388 " (.part file). It toggles between color-filled triangles and colored\n");
389 printf(
390 " edges (as the -f switch does). Filled triangles look better when all\n");
391 printf(
392 " triangles are reasonably large; colored edges look better when there\n");
393 printf(" are very small triangles present.\n");
394 printf(
395 " PS: Creates a PostScript file containing the image you are viewing. If\n"
396 );
397 printf(
398 " the -b switch is selected, all PostScript output will be black and\n");
399 printf(
400 " white; otherwise, .part.ps and .adj.ps files will be color, independent\n"
401 );
402 printf(
403 " of whether you are using a color monitor. Normally the output will\n");
404 printf(
405 " preserve the properties of the image you see on the screen, including\n");
406 printf(
407 " zoom and line width; however, if black and white output is selected (-b\n"
408 );
409 printf(
410 " switch), partitioned meshes will always be drawn exploded. The output\n"
411 );
412 printf(
413 " file name depends on the image being viewed. If you want several\n");
414 printf(
415 " different snapshots (zooming in on different parts) of the same object,\n"
416 );
417 printf(
418 " you'll have to rename each file after Show Me creates it so that it\n");
419 printf(" isn't overwritten by the next snapshot.\n");
420 printf(
421 " EPS: Creates an encapsulated PostScript file, suitable for inclusion in\n"
422 );
423 printf(
424 " documents. Otherwise, this button is just like the PS button. (The\n");
425 printf(
426 " main difference is that .eps files lack a `showpage' command at the\n");
427 printf(" end.)\n\n");
428 printf(
429 " There are two nearly-identical rows of buttons that load different images\n"
430 );
431 printf(" from disk. Each row contains the following buttons:\n\n");
432 printf(" node: Loads a .node file.\n");
433 printf(
434 " poly: Loads a .poly file (and possibly an associated .node file).\n");
435 printf(" ele: Loads an .ele file (and associated .node file).\n");
436 printf(" edge: Loads an .edge file (and associated .node file).\n");
437 printf(
438 " part: Loads a .part file (and associated .node and .ele files).\n");
439 printf(
440 " adj: Loads an .adj file (and associated .node, .ele, and .part files).\n");
441 printf(" voro: Loads a .v.node and .v.edge file for a Voronoi diagram.\n");
442 printf("\n");
443 printf(
444 " Each row represents a different iteration number of the geometry files.\n");
445 printf(
446 " For a full explanation of iteration numbers, read the instructions for\n");
447 printf(
448 " Triangle. Briefly, iteration numbers are used to allow a user to easily\n"
449 );
450 printf(
451 " represent a sequence of related triangulations. Iteration numbers are\n");
452 printf(
453 " used in the names of geometry files; for instance, mymesh.3.ele is a\n");
454 printf(
455 " triangle file with iteration number three, and mymesh.ele has an implicit\n"
456 );
457 printf(" iteration number of zero.\n\n");
458 printf(
459 " The control buttons at the right end of each row display the two\n");
460 printf(
461 " iterations currently under view. These buttons can be clicked to\n");
462 printf(
463 " increase or decrease the iteration numbers, and thus conveniently view\n");
464 printf(" a sequence of meshes.\n\n");
465 printf(
466 " Show Me keeps each file in memory after loading it, but you can force\n");
467 printf(
468 " Show Me to reread a set of files (for one iteration number) by reclicking\n"
469 );
470 printf(
471 " the button that corresponds to the current image. This is convenient if\n"
472 );
473 printf(" you have changed a geometry file.\n\n");
474 printf("File Formats:\n\n");
475 printf(
476 " All files may contain comments prefixed by the character '#'. Points,\n");
477 printf(
478 " segments, holes, triangles, edges, and subdomains must be numbered\n");
479 printf(
480 " consecutively, starting from either 1 or 0. Whichever you choose, all\n");
481 printf(
482 " input files must be consistent (for any single iteration number); if the\n"
483 );
484 printf(
485 " nodes are numbered from 1, so must be all other objects. Show Me\n");
486 printf(
487 " automatically detects your choice while reading a .node (or .poly) file.\n"
488 );
489 printf(" Examples of these file formats are given below.\n\n");
490 printf(" .node files:\n");
491 printf(
492 " First line: <# of points> <dimension (must be 2)> <# of attributes>\n");
493 printf(
494 " <# of boundary markers (0 or 1)>\n"
495 );
496 printf(
497 " Remaining lines: <point #> <x> <y> [attributes] [boundary marker]\n");
498 printf("\n");
499 printf(
500 " The attributes, which are typically floating-point values of physical\n");
501 printf(
502 " quantities (such as mass or conductivity) associated with the nodes of\n"
503 );
504 printf(
505 " a finite element mesh, are ignored by Show Me. Show Me also ignores\n");
506 printf(
507 " boundary markers. See the instructions for Triangle to find out what\n");
508 printf(" attributes and boundary markers are.\n\n");
509 printf(" .poly files:\n");
510 printf(
511 " First line: <# of points> <dimension (must be 2)> <# of attributes>\n");
512 printf(
513 " <# of boundary markers (0 or 1)>\n"
514 );
515 printf(
516 " Following lines: <point #> <x> <y> [attributes] [boundary marker]\n");
517 printf(" One line: <# of segments> <# of boundary markers (0 or 1)>\n");
518 printf(
519 " Following lines: <segment #> <endpoint> <endpoint> [boundary marker]\n");
520 printf(" One line: <# of holes>\n");
521 printf(" Following lines: <hole #> <x> <y>\n");
522 printf(" [Optional additional lines that are ignored]\n\n");
523 printf(
524 " A .poly file represents a Planar Straight Line Graph (PSLG), an idea\n");
525 printf(
526 " familiar to computational geometers. By definition, a PSLG is just a\n");
527 printf(
528 " list of points and edges. A .poly file also contains some additional\n");
529 printf(" information.\n\n");
530 printf(
531 " The first section lists all the points, and is identical to the format\n"
532 );
533 printf(
534 " of .node files. <# of points> may be set to zero to indicate that the\n"
535 );
536 printf(
537 " points are listed in a separate .node file; .poly files produced by\n");
538 printf(
539 " Triangle always have this format. When Show Me reads such a file, it\n");
540 printf(" also reads the corresponding .node file.\n\n");
541 printf(
542 " The second section lists the segments. Segments are edges whose\n");
543 printf(
544 " presence in a triangulation produced from the PSLG is enforced. Each\n");
545 printf(
546 " segment is specified by listing the indices of its two endpoints. This\n"
547 );
548 printf(
549 " means that its endpoints must be included in the point list. Each\n");
550 printf(
551 " segment, like each point, may have a boundary marker, which is ignored\n"
552 );
553 printf(" by Show Me.\n\n");
554 printf(
555 " The third section lists holes and concavities that are desired in any\n");
556 printf(
557 " triangulation generated from the PSLG. Holes are specified by\n");
558 printf(" identifying a point inside each hole.\n\n");
559 printf(" .ele files:\n");
560 printf(
561 " First line: <# of triangles> <points per triangle> <# of attributes>\n");
562 printf(
563 " Remaining lines: <triangle #> <point> <point> <point> ... [attributes]\n"
564 );
565 printf("\n");
566 printf(
567 " Points are indices into the corresponding .node file. Show Me ignores\n"
568 );
569 printf(
570 " all but the first three points of each triangle; these should be the\n");
571 printf(
572 " corners listed in counterclockwise order around the triangle. The\n");
573 printf(" attributes are ignored by Show Me.\n\n");
574 printf(" .edge files:\n");
575 printf(" First line: <# of edges> <# of boundary markers (0 or 1)>\n");
576 printf(
577 " Following lines: <edge #> <endpoint> <endpoint> [boundary marker]\n");
578 printf("\n");
579 printf(
580 " Endpoints are indices into the corresponding .node file. The boundary\n"
581 );
582 printf(" markers are ignored by Show Me.\n\n");
583 printf(
584 " In Voronoi diagrams, one also finds a special kind of edge that is an\n");
585 printf(
586 " infinite ray with only one endpoint. For these edges, a different\n");
587 printf(" format is used:\n\n");
588 printf(" <edge #> <endpoint> -1 <direction x> <direction y>\n\n");
589 printf(
590 " The `direction' is a floating-point vector that indicates the direction\n"
591 );
592 printf(" of the infinite ray.\n\n");
593 printf(" .part files:\n");
594 printf(" First line: <# of triangles> <# of subdomains>\n");
595 printf(" Remaining lines: <triangle #> <subdomain #>\n\n");
596 printf(
597 " The set of triangles is partitioned by a .part file; each triangle is\n");
598 printf(" mapped to a subdomain.\n\n");
599 printf(" .adj files:\n");
600 printf(" First line: <# of subdomains>\n");
601 printf(" Remaining lines: <adjacency matrix entry>\n\n");
602 printf(
603 " An .adj file represents adjacencies between subdomains (presumably\n");
604 printf(" computed by a partitioner). The first line is followed by\n");
605 printf(
606 " (subdomains X subdomains) lines, each containing one entry of the\n");
607 printf(
608 " adjacency matrix. A nonzero entry indicates that two subdomains are\n");
609 printf(" adjacent (share a point).\n\n");
610 printf("Example:\n\n");
611 printf(
612 " Here is a sample file `box.poly' describing a square with a square hole:\n"
613 );
614 printf("\n");
615 printf(
616 " # A box with eight points in 2D, no attributes, no boundary marker.\n");
617 printf(" 8 2 0 0\n");
618 printf(" # Outer box has these vertices:\n");
619 printf(" 1 0 0\n");
620 printf(" 2 0 3\n");
621 printf(" 3 3 0\n");
622 printf(" 4 3 3\n");
623 printf(" # Inner square has these vertices:\n");
624 printf(" 5 1 1\n");
625 printf(" 6 1 2\n");
626 printf(" 7 2 1\n");
627 printf(" 8 2 2\n");
628 printf(" # Five segments without boundary markers.\n");
629 printf(" 5 0\n");
630 printf(" 1 1 2 # Left side of outer box.\n");
631 printf(" 2 5 7 # Segments 2 through 5 enclose the hole.\n");
632 printf(" 3 7 8\n");
633 printf(" 4 8 6\n");
634 printf(" 5 6 5\n");
635 printf(" # One hole in the middle of the inner square.\n");
636 printf(" 1\n");
637 printf(" 1 1.5 1.5\n\n");
638 printf(
639 " After this PSLG is triangulated by Triangle, the resulting triangulation\n"
640 );
641 printf(
642 " consists of a .node and .ele file. Here is the former, `box.1.node',\n");
643 printf(" which duplicates the points of the PSLG:\n\n");
644 printf(" 8 2 0 0\n");
645 printf(" 1 0 0\n");
646 printf(" 2 0 3\n");
647 printf(" 3 3 0\n");
648 printf(" 4 3 3\n");
649 printf(" 5 1 1\n");
650 printf(" 6 1 2\n");
651 printf(" 7 2 1\n");
652 printf(" 8 2 2\n");
653 printf(" # Generated by triangle -pcBev box\n");
654 printf("\n");
655 printf(" Here is the triangulation file, `box.1.ele'.\n");
656 printf("\n");
657 printf(" 8 3 0\n");
658 printf(" 1 1 5 6\n");
659 printf(" 2 5 1 3\n");
660 printf(" 3 2 6 8\n");
661 printf(" 4 6 2 1\n");
662 printf(" 5 7 3 4\n");
663 printf(" 6 3 7 5\n");
664 printf(" 7 8 4 2\n");
665 printf(" 8 4 8 7\n");
666 printf(" # Generated by triangle -pcBev box\n\n");
667 printf(" Here is the edge file for the triangulation, `box.1.edge'.\n\n");
668 printf(" 16 0\n");
669 printf(" 1 1 5\n");
670 printf(" 2 5 6\n");
671 printf(" 3 6 1\n");
672 printf(" 4 1 3\n");
673 printf(" 5 3 5\n");
674 printf(" 6 2 6\n");
675 printf(" 7 6 8\n");
676 printf(" 8 8 2\n");
677 printf(" 9 2 1\n");
678 printf(" 10 7 3\n");
679 printf(" 11 3 4\n");
680 printf(" 12 4 7\n");
681 printf(" 13 7 5\n");
682 printf(" 14 8 4\n");
683 printf(" 15 4 2\n");
684 printf(" 16 8 7\n");
685 printf(" # Generated by triangle -pcBev box\n");
686 printf("\n");
687 printf(
688 " Here's a file `box.1.part' that partitions the mesh into four subdomains.\n"
689 );
690 printf("\n");
691 printf(" 8 4\n");
692 printf(" 1 3\n");
693 printf(" 2 3\n");
694 printf(" 3 4\n");
695 printf(" 4 4\n");
696 printf(" 5 1\n");
697 printf(" 6 1\n");
698 printf(" 7 2\n");
699 printf(" 8 2\n");
700 printf(" # Generated by slice -s4 box.1\n\n");
701 printf(
702 " Here's a file `box.1.adj' that represents the resulting adjacencies.\n");
703 printf("\n");
704 printf(" 4\n");
705 printf(" 9\n");
706 printf(" 2\n");
707 printf(" 2\n");
708 printf(" 0\n");
709 printf(" 2\n");
710 printf(" 9\n");
711 printf(" 0\n");
712 printf(" 2\n");
713 printf(" 2\n");
714 printf(" 0\n");
715 printf(" 9\n");
716 printf(" 2\n");
717 printf(" 0\n");
718 printf(" 2\n");
719 printf(" 2\n");
720 printf(" 9\n");
721 printf("\n");
722 printf("Display Speed:\n");
723 printf("\n");
724 printf(
725 " It is worthwhile to note that .edge files typically plot and print twice\n"
726 );
727 printf(
728 " as quickly as .ele files, because .ele files cause each internal edge to\n"
729 );
730 printf(
731 " be drawn twice. For the same reason, PostScript files created from edge\n"
732 );
733 printf(" sets are smaller than those created from triangulations.\n\n");
734 printf("Show Me on the Web:\n\n");
735 printf(
736 " To see an illustrated, updated version of these instructions, check out\n");
737 printf("\n");
738 printf(" http://www.cs.cmu.edu/~quake/showme.html\n");
739 printf("\n");
740 printf("A Brief Plea:\n");
741 printf("\n");
742 printf(
743 " If you use Show Me (or Triangle), and especially if you use it to\n");
744 printf(
745 " accomplish real work, I would like very much to hear from you. A short\n");
746 printf(
747 " letter or email (to jrs@cs.cmu.edu) describing how you use Show Me (and\n");
748 printf(
749 " its sister programs) will mean a lot to me. The more people I know\n");
750 printf(
751 " are using my programs, the more easily I can justify spending time on\n");
752 printf(
753 " improvements, which in turn will benefit you. Also, I can put you\n");
754 printf(
755 " on a list to receive email whenever new versions are available.\n");
756 printf("\n");
757 printf(
758 " If you use a PostScript file generated by Show Me in a publication,\n");
759 printf(" please include an acknowledgment as well.\n\n");
760 exit(0);
761 }
762
set_filenames(filename,lowermeshnumber)763 void set_filenames(filename, lowermeshnumber)
764 char *filename;
765 int lowermeshnumber;
766 {
767 char numberstring[100];
768 int i;
769
770 for (i = 0; i < 2; i++) {
771 strcpy(nodefilename[i], filename);
772 strcpy(polyfilename[i], filename);
773 strcpy(elefilename[i], filename);
774 strcpy(edgefilename[i], filename);
775 strcpy(partfilename[i], filename);
776 strcpy(adjfilename[i], filename);
777 strcpy(vnodefilename[i], filename);
778 strcpy(vedgefilename[i], filename);
779
780 if (lowermeshnumber + i > 0) {
781 sprintf(numberstring, ".%d", lowermeshnumber + i);
782 strcat(nodefilename[i], numberstring);
783 strcat(polyfilename[i], numberstring);
784 strcat(elefilename[i], numberstring);
785 strcat(edgefilename[i], numberstring);
786 strcat(partfilename[i], numberstring);
787 strcat(adjfilename[i], numberstring);
788 strcat(vnodefilename[i], numberstring);
789 strcat(vedgefilename[i], numberstring);
790 }
791
792 strcat(nodefilename[i], ".node");
793 strcat(polyfilename[i], ".poly");
794 strcat(elefilename[i], ".ele");
795 strcat(edgefilename[i], ".edge");
796 strcat(partfilename[i], ".part");
797 strcat(adjfilename[i], ".adj");
798 strcat(vnodefilename[i], ".v.node");
799 strcat(vedgefilename[i], ".v.edge");
800 }
801 }
802
803 #if 1 /* This function is already in netlib.lib, see triangle.c */
804 void parsecommandline(int argc, char **argv);
805 #else
parsecommandline(argc,argv)806 void parsecommandline(argc, argv)
807 int argc;
808 char **argv;
809 {
810 int increment;
811 int meshnumber;
812 int i, j;
813
814 quiet = 0;
815 fillelem = 0;
816 line_width = 1;
817 bw_ps = 0;
818 start_image = ELE;
819 filename[0] = '\0';
820 for (i = 1; i < argc; i++) {
821 if (argv[i][0] == '-') {
822 for (j = 1; argv[i][j] != '\0'; j++) {
823 if (argv[i][j] == 'f') {
824 fillelem = 1;
825 }
826 if (argv[i][j] == 'w') {
827 if ((argv[i][j + 1] >= '1') && (argv[i][j + 1] <= '9')) {
828 line_width = 0;
829 while ((argv[i][j + 1] >= '0') && (argv[i][j + 1] <= '9')) {
830 j++;
831 line_width = line_width * 10 + (int) (argv[i][j] - '0');
832 }
833 if (line_width > 100) {
834 printf("Error: Line width cannot exceed 100.\n");
835 line_width = 1;
836 }
837 }
838 }
839 if (argv[i][j] == 'b') {
840 bw_ps = 1;
841 }
842 if (argv[i][j] == 'Q') {
843 quiet = 1;
844 }
845 if ((argv[i][j] == 'h') || (argv[i][j] == 'H') ||
846 (argv[i][j] == '?')) {
847 info();
848 }
849 }
850 } else {
851 strcpy(filename, argv[i]);
852 }
853 }
854 if (filename[0] == '\0') {
855 syntax();
856 }
857 if (!strcmp(&filename[strlen(filename) - 5], ".node")) {
858 filename[strlen(filename) - 5] = '\0';
859 start_image = NODE;
860 }
861 if (!strcmp(&filename[strlen(filename) - 5], ".poly")) {
862 filename[strlen(filename) - 5] = '\0';
863 start_image = POLY;
864 }
865 if (!strcmp(&filename[strlen(filename) - 4], ".ele")) {
866 filename[strlen(filename) - 4] = '\0';
867 start_image = ELE;
868 }
869 if (!strcmp(&filename[strlen(filename) - 5], ".edge")) {
870 filename[strlen(filename) - 5] = '\0';
871 start_image = EDGE;
872 }
873 if (!strcmp(&filename[strlen(filename) - 5], ".part")) {
874 filename[strlen(filename) - 5] = '\0';
875 start_image = PART;
876 }
877 if (!strcmp(&filename[strlen(filename) - 4], ".adj")) {
878 filename[strlen(filename) - 4] = '\0';
879 start_image = ADJ;
880 }
881
882 increment = 0;
883 j = 1;
884 while (filename[j] != '\0') {
885 if ((filename[j] == '.') && (filename[j + 1] != '\0')) {
886 increment = j + 1;
887 }
888 j++;
889 }
890 meshnumber = 0;
891 if (increment > 0) {
892 j = increment;
893 do {
894 if ((filename[j] >= '0') && (filename[j] <= '9')) {
895 meshnumber = meshnumber * 10 + (int) (filename[j] - '0');
896 } else {
897 increment = 0;
898 }
899 j++;
900 } while (filename[j] != '\0');
901 }
902 if (increment > 0) {
903 filename[increment - 1] = '\0';
904 }
905
906 if (meshnumber == 0) {
907 start_inc = 0;
908 loweriteration = 0;
909 } else {
910 start_inc = 1;
911 loweriteration = meshnumber - 1;
912 }
913 set_filenames(filename, loweriteration);
914 }
915 #endif /* 0 */
916
free_inc(inc)917 void free_inc(inc)
918 int inc;
919 {
920 if (loaded[inc][NODE]) {
921 free(nodeptr[inc]);
922 }
923 if (loaded[inc][POLY]) {
924 if (polynodes[inc] > 0) {
925 free(polynodeptr[inc]);
926 }
927 free(polyedgeptr[inc]);
928 free(polyholeptr[inc]);
929 }
930 if (loaded[inc][ELE]) {
931 free(eleptr[inc]);
932 }
933 if (loaded[inc][PART]) {
934 free(partpart[inc]);
935 free(partcenter[inc]);
936 free(partshift[inc]);
937 }
938 if (loaded[inc][EDGE]) {
939 free(edgeptr[inc]);
940 free(normptr[inc]);
941 }
942 if (loaded[inc][ADJ]) {
943 free(adjptr[inc]);
944 }
945 if (loaded[inc][VORO]) {
946 free(vnodeptr[inc]);
947 free(vedgeptr[inc]);
948 free(vnormptr[inc]);
949 }
950 }
951
move_inc(inc)952 void move_inc(inc)
953 int inc;
954 {
955 int i;
956
957 free_inc(1 - inc);
958 for (i = 0; i < IMAGE_TYPES; i++) {
959 loaded[1 - inc][i] = loaded[inc][i];
960 loaded[inc][i] = 0;
961 xlo[1 - inc][i] = xlo[inc][i];
962 ylo[1 - inc][i] = ylo[inc][i];
963 xhi[1 - inc][i] = xhi[inc][i];
964 yhi[1 - inc][i] = yhi[inc][i];
965 }
966 nodes[1 - inc] = nodes[inc];
967 node_dim[1 - inc] = node_dim[inc];
968 nodeptr[1 - inc] = nodeptr[inc];
969 polynodes[1 - inc] = polynodes[inc];
970 poly_dim[1 - inc] = poly_dim[inc];
971 polyedges[1 - inc] = polyedges[inc];
972 polyholes[1 - inc] = polyholes[inc];
973 polynodeptr[1 - inc] = polynodeptr[inc];
974 polyedgeptr[1 - inc] = polyedgeptr[inc];
975 polyholeptr[1 - inc] = polyholeptr[inc];
976 elems[1 - inc] = elems[inc];
977 ele_corners[1 - inc] = ele_corners[inc];
978 eleptr[1 - inc] = eleptr[inc];
979 edges[1 - inc] = edges[inc];
980 edgeptr[1 - inc] = edgeptr[inc];
981 normptr[1 - inc] = normptr[inc];
982 subdomains[1 - inc] = subdomains[inc];
983 partpart[1 - inc] = partpart[inc];
984 partcenter[1 - inc] = partcenter[inc];
985 partshift[1 - inc] = partshift[inc];
986 adjsubdomains[1 - inc] = adjsubdomains[inc];
987 adjptr[1 - inc] = adjptr[inc];
988 vnodes[1 - inc] = vnodes[inc];
989 vnode_dim[1 - inc] = vnode_dim[inc];
990 vnodeptr[1 - inc] = vnodeptr[inc];
991 vedges[1 - inc] = vedges[inc];
992 vedgeptr[1 - inc] = vedgeptr[inc];
993 vnormptr[1 - inc] = vnormptr[inc];
994 firstnumber[1 - inc] = firstnumber[inc];
995 firstnumber[inc] = -1;
996 }
997
unload_inc(inc)998 void unload_inc(inc)
999 int inc;
1000 {
1001 int i;
1002
1003 current_image = NOTHING;
1004 for (i = 0; i < IMAGE_TYPES; i++) {
1005 loaded[inc][i] = 0;
1006 firstnumber[inc] = -1;
1007 }
1008 }
1009
showme_init()1010 void showme_init()
1011 {
1012 current_image = NOTHING;
1013 current_inc = 0;
1014 explosion = STARTEXPLOSION;
1015 unload_inc(0);
1016 unload_inc(1);
1017 }
1018
readline(string,infile,infilename)1019 char *readline(string, infile, infilename)
1020 char *string;
1021 FILE *infile;
1022 char *infilename;
1023 {
1024 char *result;
1025
1026 do {
1027 result = fgets(string, INPUTLINESIZE, infile);
1028 if (result == (char *) NULL) {
1029 printf(" Error: Unexpected end of file in %s.\n",
1030 infilename);
1031 exit(1);
1032 }
1033 while ((*result != '\0') && (*result != '#')
1034 && (*result != '.') && (*result != '+') && (*result != '-')
1035 && ((*result < '0') || (*result > '9'))) {
1036 result++;
1037 }
1038 } while ((*result == '#') || (*result == '\0'));
1039 return result;
1040 }
1041
findfield(string)1042 char *findfield(string)
1043 char *string;
1044 {
1045 char *result;
1046
1047 result = string;
1048 while ((*result != '\0') && (*result != '#')
1049 && (*result != ' ') && (*result != '\t')) {
1050 result++;
1051 }
1052 while ((*result != '\0') && (*result != '#')
1053 && (*result != '.') && (*result != '+') && (*result != '-')
1054 && ((*result < '0') || (*result > '9'))) {
1055 result++;
1056 }
1057 if (*result == '#') {
1058 *result = '\0';
1059 }
1060 return result;
1061 }
1062
load_node(fname,firstnumber,nodes,dim,ptr,xmin,ymin,xmax,ymax)1063 int load_node(fname, firstnumber, nodes, dim, ptr, xmin, ymin, xmax, ymax)
1064 char *fname;
1065 int *firstnumber;
1066 int *nodes;
1067 int *dim;
1068 REAL **ptr;
1069 REAL *xmin;
1070 REAL *ymin;
1071 REAL *xmax;
1072 REAL *ymax;
1073 {
1074 FILE *infile;
1075 char inputline[INPUTLINESIZE];
1076 char *stringptr;
1077 int extras;
1078 int nodemarks;
1079 int index;
1080 int nodenumber;
1081 int i, j;
1082 int smallerr;
1083 REAL x, y;
1084
1085 *xmin = *ymin = 0.0;
1086 *xmax = *ymax = 1.0;
1087 if (!quiet) {
1088 printf("Opening %s.\n", fname);
1089 }
1090 infile = fopen(fname, "r");
1091 if (infile == (FILE *) NULL) {
1092 printf(" Error: Cannot access file %s.\n", fname);
1093 return 1;
1094 }
1095 stringptr = readline(inputline, infile, fname);
1096 *nodes = (int) strtol (stringptr, &stringptr, 0);
1097 if (*nodes < 3) {
1098 printf(" Error: %s contains %d points.\n", fname, *nodes);
1099 return 1;
1100 }
1101 stringptr = findfield(stringptr);
1102 if (*stringptr == '\0') {
1103 *dim = 2;
1104 } else {
1105 *dim = (int) strtol (stringptr, &stringptr, 0);
1106 }
1107 if (*dim < 1) {
1108 printf(" Error: %s has dimensionality %d.\n", fname, *dim);
1109 return 1;
1110 }
1111 if (*dim != 2) {
1112 printf(" I only understand two-dimensional meshes.\n");
1113 return 1;
1114 }
1115 stringptr = findfield(stringptr);
1116 if (*stringptr == '\0') {
1117 extras = 0;
1118 } else {
1119 extras = (int) strtol (stringptr, &stringptr, 0);
1120 }
1121 if (extras < 0) {
1122 printf(" Error: %s has negative value for number of attributes.\n",
1123 fname);
1124 return 1;
1125 }
1126 stringptr = findfield(stringptr);
1127 if (*stringptr == '\0') {
1128 nodemarks = 0;
1129 } else {
1130 nodemarks = (int) strtol (stringptr, &stringptr, 0);
1131 }
1132 if (nodemarks < 0) {
1133 printf(" Warning: %s has negative value for number of point markers.\n",
1134 fname);
1135 }
1136 if (nodemarks > 1) {
1137 printf(
1138 " Warning: %s has value greater than one for number of point markers.\n",
1139 fname);
1140 }
1141 *ptr = (REAL *) malloc((*nodes + 1) * *dim * sizeof(REAL));
1142 if (*ptr == (REAL *) NULL) {
1143 printf(" Out of memory.\n");
1144 return 1;
1145 }
1146 index = *dim;
1147 smallerr = 1;
1148 for (i = 0; i < *nodes; i++) {
1149 stringptr = readline(inputline, infile, fname);
1150 nodenumber = (int) strtol (stringptr, &stringptr, 0);
1151 if ((i == 0) && (*firstnumber == -1)) {
1152 if (nodenumber == 0) {
1153 *firstnumber = 0;
1154 } else {
1155 *firstnumber = 1;
1156 }
1157 }
1158 if ((nodenumber != *firstnumber + i) && (smallerr)) {
1159 printf(" Warning: Points in %s are not numbered correctly\n", fname);
1160 printf(" (starting with point %d).\n", *firstnumber + i);
1161 smallerr = 0;
1162 }
1163 for (j = 0; j < *dim; j++) {
1164 stringptr = findfield(stringptr);
1165 if (*stringptr == '\0') {
1166 printf("Error: Point %d is missing a coordinate in %s.\n",
1167 *firstnumber + i, fname);
1168 free(*ptr);
1169 return 1;
1170 }
1171 (*ptr)[index++] = (REAL) strtod(stringptr, &stringptr);
1172 }
1173 }
1174 fclose(infile);
1175 index = *dim;
1176 *xmin = *xmax = (*ptr)[index];
1177 *ymin = *ymax = (*ptr)[index + 1];
1178 for (i = 2; i <= *nodes; i++) {
1179 index += *dim;
1180 x = (*ptr)[index];
1181 y = (*ptr)[index + 1];
1182 if (x < *xmin) {
1183 *xmin = x;
1184 }
1185 if (y < *ymin) {
1186 *ymin = y;
1187 }
1188 if (x > *xmax) {
1189 *xmax = x;
1190 }
1191 if (y > *ymax) {
1192 *ymax = y;
1193 }
1194 }
1195 if (*xmin == *xmax) {
1196 *xmin -= 0.5;
1197 *xmax += 0.5;
1198 }
1199 if (*ymin == *ymax) {
1200 *ymin -= 0.5;
1201 *ymax += 0.5;
1202 }
1203 return 0;
1204 }
1205
load_poly(inc,fname,firstnumber,pnodes,dim,edges,holes,nodeptr,edgeptr,holeptr,xmin,ymin,xmax,ymax)1206 int load_poly(inc, fname, firstnumber, pnodes, dim, edges, holes, nodeptr,
1207 edgeptr, holeptr, xmin, ymin, xmax, ymax)
1208 int inc;
1209 char *fname;
1210 int *firstnumber;
1211 int *pnodes;
1212 int *dim;
1213 int *edges;
1214 int *holes;
1215 REAL **nodeptr;
1216 int **edgeptr;
1217 REAL **holeptr;
1218 REAL *xmin;
1219 REAL *ymin;
1220 REAL *xmax;
1221 REAL *ymax;
1222 {
1223 FILE *infile;
1224 char inputline[INPUTLINESIZE];
1225 char *stringptr;
1226 int extras;
1227 int nodemarks;
1228 int segmentmarks;
1229 int index;
1230 int nodenumber, edgenumber, holenumber;
1231 int maxnode;
1232 int i, j;
1233 int smallerr;
1234 REAL x, y;
1235
1236 if (!quiet) {
1237 printf("Opening %s.\n", fname);
1238 }
1239 infile = fopen(fname, "r");
1240 if (infile == (FILE *) NULL) {
1241 printf(" Error: Cannot access file %s.\n", fname);
1242 return 1;
1243 }
1244 stringptr = readline(inputline, infile, fname);
1245 *pnodes = (int) strtol (stringptr, &stringptr, 0);
1246 if (*pnodes == 0) {
1247 if (!loaded[inc][NODE]) {
1248 if (load_image(inc, NODE)) {
1249 return 1;
1250 }
1251 }
1252 maxnode = nodes[inc];
1253 *xmin = xlo[inc][NODE];
1254 *ymin = ylo[inc][NODE];
1255 *xmax = xhi[inc][NODE];
1256 *ymax = yhi[inc][NODE];
1257 } else {
1258 if (*pnodes < 1) {
1259 printf(" Error: %s contains %d points.\n", fname, *pnodes);
1260 return 1;
1261 }
1262 maxnode = *pnodes;
1263 }
1264 stringptr = findfield(stringptr);
1265 if (*stringptr == '\0') {
1266 *dim = 2;
1267 } else {
1268 *dim = (int) strtol (stringptr, &stringptr, 0);
1269 }
1270 if (*dim < 1) {
1271 printf(" Error: %s has dimensionality %d.\n", fname, *dim);
1272 return 1;
1273 }
1274 if (*dim != 2) {
1275 printf(" I only understand two-dimensional meshes.\n");
1276 return 1;
1277 }
1278 stringptr = findfield(stringptr);
1279 if (*stringptr == '\0') {
1280 extras = 0;
1281 } else {
1282 extras = (int) strtol (stringptr, &stringptr, 0);
1283 }
1284 if (extras < 0) {
1285 printf(" Error: %s has negative value for number of attributes.\n",
1286 fname);
1287 return 1;
1288 }
1289 stringptr = findfield(stringptr);
1290 if (*stringptr == '\0') {
1291 nodemarks = 0;
1292 } else {
1293 nodemarks = (int) strtol (stringptr, &stringptr, 0);
1294 }
1295 if (nodemarks < 0) {
1296 printf(" Warning: %s has negative value for number of point markers.\n",
1297 fname);
1298 }
1299 if (nodemarks > 1) {
1300 printf(
1301 " Warning: %s has value greater than one for number of point markers.\n",
1302 fname);
1303 }
1304 if (*pnodes > 0) {
1305 *nodeptr = (REAL *) malloc((*pnodes + 1) * *dim * sizeof(REAL));
1306 if (*nodeptr == (REAL *) NULL) {
1307 printf(" Out of memory.\n");
1308 return 1;
1309 }
1310 index = *dim;
1311 smallerr = 1;
1312 for (i = 0; i < *pnodes; i++) {
1313 stringptr = readline(inputline, infile, fname);
1314 nodenumber = (int) strtol (stringptr, &stringptr, 0);
1315 if ((i == 0) && (*firstnumber == -1)) {
1316 if (nodenumber == 0) {
1317 *firstnumber = 0;
1318 } else {
1319 *firstnumber = 1;
1320 }
1321 }
1322 if ((nodenumber != *firstnumber + i) && (smallerr)) {
1323 printf(" Warning: Points in %s are not numbered correctly.\n",
1324 fname);
1325 printf(" (starting with point %d).\n", *firstnumber + i);
1326 smallerr = 0;
1327 }
1328 for (j = 0; j < *dim; j++) {
1329 stringptr = findfield(stringptr);
1330 if (*stringptr == '\0') {
1331 printf("Error: Point %d is missing a coordinate in %s.\n",
1332 *firstnumber + i, fname);
1333 free(*nodeptr);
1334 return 1;
1335 }
1336 (*nodeptr)[index++] = (REAL) strtod(stringptr, &stringptr);
1337 }
1338 }
1339 }
1340 stringptr = readline(inputline, infile, fname);
1341 *edges = (int) strtol (stringptr, &stringptr, 0);
1342 if (*edges < 0) {
1343 printf(" Error: %s contains %d segments.\n", fname, *edges);
1344 free(*nodeptr);
1345 return 1;
1346 }
1347 stringptr = findfield(stringptr);
1348 if (*stringptr == '\0') {
1349 segmentmarks = 0;
1350 } else {
1351 segmentmarks = (int) strtol (stringptr, &stringptr, 0);
1352 }
1353 if (segmentmarks < 0) {
1354 printf(" Error: %s has negative value for number of segment markers.\n",
1355 fname);
1356 free(*nodeptr);
1357 return 1;
1358 }
1359 if (segmentmarks > 1) {
1360 printf(
1361 " Error: %s has value greater than one for number of segment markers.\n",
1362 fname);
1363 free(*nodeptr);
1364 return 1;
1365 }
1366 *edgeptr = (int *) malloc(((*edges + 1) << 1) * sizeof(int));
1367 if (*edgeptr == (int *) NULL) {
1368 printf(" Out of memory.\n");
1369 free(*nodeptr);
1370 return 1;
1371 }
1372 index = 2;
1373 smallerr = 1;
1374 for (i = *firstnumber; i < *firstnumber + *edges; i++) {
1375 stringptr = readline(inputline, infile, fname);
1376 edgenumber = (int) strtol (stringptr, &stringptr, 0);
1377 if ((edgenumber != i) && (smallerr)) {
1378 printf(" Warning: Segments in %s are not numbered correctly.\n",
1379 fname);
1380 printf(" (starting with segment %d).\n", i);
1381 smallerr = 0;
1382 }
1383 stringptr = findfield(stringptr);
1384 if (*stringptr == '\0') {
1385 printf("Error: Segment %d is missing its endpoints in %s.\n", i, fname);
1386 free(*nodeptr);
1387 free(*edgeptr);
1388 return 1;
1389 }
1390 (*edgeptr)[index] = (int) strtol (stringptr, &stringptr, 0) + 1 -
1391 *firstnumber;
1392 if (((*edgeptr)[index] < 1) || ((*edgeptr)[index] > maxnode)) {
1393 printf("Error: Segment %d has invalid endpoint in %s.\n", i, fname);
1394 return 1;
1395 }
1396 stringptr = findfield(stringptr);
1397 if (*stringptr == '\0') {
1398 printf("Error: Segment %d is missing an endpoint in %s.\n", i, fname);
1399 free(*nodeptr);
1400 free(*edgeptr);
1401 return 1;
1402 }
1403 (*edgeptr)[index + 1] = (int) strtol (stringptr, &stringptr, 0) + 1 -
1404 *firstnumber;
1405 if (((*edgeptr)[index + 1] < 1) || ((*edgeptr)[index + 1] > maxnode)) {
1406 printf("Error: Segment %d has invalid endpoint in %s.\n", i, fname);
1407 return 1;
1408 }
1409 index += 2;
1410 }
1411 stringptr = readline(inputline, infile, fname);
1412 *holes = (int) strtol (stringptr, &stringptr, 0);
1413 if (*holes < 0) {
1414 printf(" Error: %s contains %d holes.\n", fname, *holes);
1415 free(*nodeptr);
1416 free(*edgeptr);
1417 return 1;
1418 }
1419 *holeptr = (REAL *) malloc((*holes + 1) * *dim * sizeof(REAL));
1420 if (*holeptr == (REAL *) NULL) {
1421 printf(" Out of memory.\n");
1422 free(*nodeptr);
1423 free(*edgeptr);
1424 return 1;
1425 }
1426 index = *dim;
1427 smallerr = 1;
1428 for (i = *firstnumber; i < *firstnumber + *holes; i++) {
1429 stringptr = readline(inputline, infile, fname);
1430 holenumber = (int) strtol (stringptr, &stringptr, 0);
1431 if ((holenumber != i) && (smallerr)) {
1432 printf(" Warning: Holes in %s are not numbered correctly.\n", fname);
1433 printf(" (starting with hole %d).\n", i);
1434 smallerr = 0;
1435 }
1436 for (j = 0; j < *dim; j++) {
1437 stringptr = findfield(stringptr);
1438 if (*stringptr == '\0') {
1439 printf("Error: Hole %d is missing a coordinate in %s.\n", i,
1440 fname);
1441 free(*nodeptr);
1442 free(*edgeptr);
1443 free(*holeptr);
1444 return 1;
1445 }
1446 (*holeptr)[index++] = (REAL) strtod(stringptr, &stringptr);
1447 }
1448 }
1449 fclose(infile);
1450 if (*pnodes > 0) {
1451 index = *dim;
1452 *xmin = *xmax = (*nodeptr)[index];
1453 *ymin = *ymax = (*nodeptr)[index + 1];
1454 for (i = 2; i <= *pnodes; i++) {
1455 index += *dim;
1456 x = (*nodeptr)[index];
1457 y = (*nodeptr)[index + 1];
1458 if (x < *xmin) {
1459 *xmin = x;
1460 }
1461 if (y < *ymin) {
1462 *ymin = y;
1463 }
1464 if (x > *xmax) {
1465 *xmax = x;
1466 }
1467 if (y > *ymax) {
1468 *ymax = y;
1469 }
1470 }
1471 }
1472 index = *dim;
1473 for (i = 1; i <= *holes; i++) {
1474 x = (*holeptr)[index];
1475 y = (*holeptr)[index + 1];
1476 if (x < *xmin) {
1477 *xmin = x;
1478 }
1479 if (y < *ymin) {
1480 *ymin = y;
1481 }
1482 if (x > *xmax) {
1483 *xmax = x;
1484 }
1485 if (y > *ymax) {
1486 *ymax = y;
1487 }
1488 index += *dim;
1489 }
1490 return 0;
1491 }
1492
load_ele(fname,firstnumber,nodes,elems,corners,ptr)1493 int load_ele(fname, firstnumber, nodes, elems, corners, ptr)
1494 char *fname;
1495 int firstnumber;
1496 int nodes;
1497 int *elems;
1498 int *corners;
1499 int **ptr;
1500 {
1501 FILE *infile;
1502 char inputline[INPUTLINESIZE];
1503 char *stringptr;
1504 int extras;
1505 int index;
1506 int elemnumber;
1507 int i, j;
1508 int smallerr;
1509
1510 if (!quiet) {
1511 printf("Opening %s.\n", fname);
1512 }
1513 infile = fopen(fname, "r");
1514 if (infile == (FILE *) NULL) {
1515 printf(" Error: Cannot access file %s.\n", fname);
1516 return 1;
1517 }
1518 stringptr = readline(inputline, infile, fname);
1519 *elems = (int) strtol (stringptr, &stringptr, 0);
1520 if (*elems < 1) {
1521 printf(" Error: %s contains %d triangles.\n", fname, *elems);
1522 return 1;
1523 }
1524 stringptr = findfield(stringptr);
1525 if (*stringptr == '\0') {
1526 *corners = 3;
1527 } else {
1528 *corners = (int) strtol (stringptr, &stringptr, 0);
1529 }
1530 if (*corners < 3) {
1531 printf(" Error: Triangles in %s have only %d corners.\n", fname,
1532 *corners);
1533 return 1;
1534 }
1535 stringptr = findfield(stringptr);
1536 if (*stringptr == '\0') {
1537 extras = 0;
1538 } else {
1539 extras = (int) strtol (stringptr, &stringptr, 0);
1540 }
1541 if (extras < 0) {
1542 printf(" Error: %s has negative value for extra fields.\n", fname);
1543 return 1;
1544 }
1545 *ptr = (int *) malloc((*elems + 1) * 3 * sizeof(int));
1546 if (*ptr == (int *) NULL) {
1547 printf(" Out of memory.\n");
1548 return 1;
1549 }
1550 index = 3;
1551 smallerr = 1;
1552 for (i = firstnumber; i < firstnumber + *elems; i++) {
1553 stringptr = readline(inputline, infile, fname);
1554 elemnumber = (int) strtol (stringptr, &stringptr, 0);
1555 if ((elemnumber != i) && (smallerr)) {
1556 printf(" Warning: Triangles in %s are not numbered correctly.\n",
1557 fname);
1558 printf(" (starting with triangle %d).\n", i);
1559 smallerr = 0;
1560 }
1561 for (j = 0; j < 3; j++) {
1562 stringptr = findfield(stringptr);
1563 if (*stringptr == '\0') {
1564 printf("Error: Triangle %d is missing a corner in %s.\n", i, fname);
1565 free(*ptr);
1566 return 1;
1567 }
1568 (*ptr)[index] = (int) strtol (stringptr, &stringptr, 0) + 1 -
1569 firstnumber;
1570 if (((*ptr)[index] < 1) || ((*ptr)[index] > nodes)) {
1571 printf("Error: Triangle %d has invalid corner in %s.\n", i, fname);
1572 return 1;
1573 }
1574 index++;
1575 }
1576 }
1577 fclose(infile);
1578 return 0;
1579 }
1580
load_edge(fname,firstnumber,nodes,edges,edgeptr,normptr)1581 int load_edge(fname, firstnumber, nodes, edges, edgeptr, normptr)
1582 char *fname;
1583 int firstnumber;
1584 int nodes;
1585 int *edges;
1586 int **edgeptr;
1587 REAL **normptr;
1588 {
1589 FILE *infile;
1590 char inputline[INPUTLINESIZE];
1591 char *stringptr;
1592 int index;
1593 int edgenumber;
1594 int edgemarks;
1595 int i;
1596 int smallerr;
1597
1598 if (!quiet) {
1599 printf("Opening %s.\n", fname);
1600 }
1601 infile = fopen(fname, "r");
1602 if (infile == (FILE *) NULL) {
1603 printf(" Error: Cannot access file %s.\n", fname);
1604 return 1;
1605 }
1606 stringptr = readline(inputline, infile, fname);
1607 *edges = (int) strtol (stringptr, &stringptr, 0);
1608 if (*edges < 1) {
1609 printf(" Error: %s contains %d edges.\n", fname, *edges);
1610 return 1;
1611 }
1612 stringptr = findfield(stringptr);
1613 if (*stringptr == '\0') {
1614 edgemarks = 0;
1615 } else {
1616 edgemarks = (int) strtol (stringptr, &stringptr, 0);
1617 }
1618 if (edgemarks < 0) {
1619 printf(" Error: %s has negative value for number of edge markers.\n",
1620 fname);
1621 return 1;
1622 }
1623 if (edgemarks > 1) {
1624 printf(
1625 " Error: %s has value greater than one for number of edge markers.\n",
1626 fname);
1627 return 1;
1628 }
1629 *edgeptr = (int *) malloc(((*edges + 1) << 1) * sizeof(int));
1630 if (*edgeptr == (int *) NULL) {
1631 printf(" Out of memory.\n");
1632 return 1;
1633 }
1634 *normptr = (REAL *) malloc(((*edges + 1) << 1) * sizeof(REAL));
1635 if (*normptr == (REAL *) NULL) {
1636 printf(" Out of memory.\n");
1637 free(*edgeptr);
1638 return 1;
1639 }
1640 index = 2;
1641 smallerr = 1;
1642 for (i = firstnumber; i < firstnumber + *edges; i++) {
1643 stringptr = readline(inputline, infile, fname);
1644 edgenumber = (int) strtol (stringptr, &stringptr, 0);
1645 if ((edgenumber != i) && (smallerr)) {
1646 printf(" Warning: Edges in %s are not numbered correctly.\n", fname);
1647 printf(" (starting with edge %d).\n", i);
1648 smallerr = 0;
1649 }
1650 stringptr = findfield(stringptr);
1651 if (*stringptr == '\0') {
1652 printf("Error: Edge %d is missing its endpoints in %s.\n", i, fname);
1653 free(*edgeptr);
1654 free(*normptr);
1655 return 1;
1656 }
1657 (*edgeptr)[index] = (int) strtol (stringptr, &stringptr, 0) + 1 -
1658 firstnumber;
1659 if (((*edgeptr)[index] < 1) || ((*edgeptr)[index] > nodes)) {
1660 printf("Error: Edge %d has invalid endpoint in %s.\n", i, fname);
1661 return 1;
1662 }
1663 stringptr = findfield(stringptr);
1664 if (*stringptr == '\0') {
1665 printf("Error: Edge %d is missing an endpoint in %s.\n", i, fname);
1666 free(*edgeptr);
1667 free(*normptr);
1668 return 1;
1669 }
1670 (*edgeptr)[index + 1] = (int) strtol (stringptr, &stringptr, 0);
1671 if ((*edgeptr)[index + 1] == -1) {
1672 stringptr = findfield(stringptr);
1673 if (*stringptr == '\0') {
1674 printf("Error: Edge %d is missing its direction in %s.\n", i, fname);
1675 free(*edgeptr);
1676 free(*normptr);
1677 return 1;
1678 }
1679 (*normptr)[index] = (REAL) strtod(stringptr, &stringptr);
1680 stringptr = findfield(stringptr);
1681 if (*stringptr == '\0') {
1682 printf("Error: Edge %d is missing a direction coordinate in %s.\n",
1683 i, fname);
1684 free(*edgeptr);
1685 free(*normptr);
1686 return 1;
1687 }
1688 (*normptr)[index + 1] = (REAL) strtod(stringptr, &stringptr);
1689 } else {
1690 (*edgeptr)[index + 1] += 1 - firstnumber;
1691 if (((*edgeptr)[index + 1] < 1) || ((*edgeptr)[index + 1] > nodes)) {
1692 printf("Error: Edge %d has invalid endpoint in %s.\n", i, fname);
1693 return 1;
1694 }
1695 }
1696 index += 2;
1697 }
1698 fclose(infile);
1699 return 0;
1700 }
1701
load_part(fname,dim,firstnumber,elems,nodeptr,eleptr,parts,partition,partcenter,partshift)1702 int load_part(fname, dim, firstnumber, elems, nodeptr, eleptr, parts,
1703 partition, partcenter, partshift)
1704 char *fname;
1705 int dim;
1706 int firstnumber;
1707 int elems;
1708 REAL *nodeptr;
1709 int *eleptr;
1710 int *parts;
1711 int **partition;
1712 REAL **partcenter;
1713 REAL **partshift;
1714 {
1715 FILE *infile;
1716 char inputline[INPUTLINESIZE];
1717 char *stringptr;
1718 int partelems;
1719 int index;
1720 int elemnumber;
1721 int i, j;
1722 int smallerr;
1723 int *partsize;
1724
1725 if (!quiet) {
1726 printf("Opening %s.\n", fname);
1727 }
1728 infile = fopen(fname, "r");
1729 if (infile == (FILE *) NULL) {
1730 printf(" Error: Cannot access file %s.\n", fname);
1731 return 1;
1732 }
1733 stringptr = readline(inputline, infile, fname);
1734 partelems = (int) strtol (stringptr, &stringptr, 0);
1735 if (partelems != elems) {
1736 printf(
1737 " Error: .ele and .part files do not agree on number of triangles.\n");
1738 return 1;
1739 }
1740 stringptr = findfield(stringptr);
1741 if (*stringptr == '\0') {
1742 *parts = 1;
1743 } else {
1744 *parts = (int) strtol (stringptr, &stringptr, 0);
1745 }
1746 if (*parts < 1) {
1747 printf(" Error: %s specifies %d subdomains.\n", fname, *parts);
1748 return 1;
1749 }
1750 *partition = (int *) malloc((elems + 1) * sizeof(int));
1751 if (*partition == (int *) NULL) {
1752 printf(" Out of memory.\n");
1753 return 1;
1754 }
1755 smallerr = 1;
1756 for (i = firstnumber; i < firstnumber + partelems; i++) {
1757 stringptr = readline(inputline, infile, fname);
1758 elemnumber = (int) strtol (stringptr, &stringptr, 0);
1759 if ((elemnumber != i) && (smallerr)) {
1760 printf(" Warning: Triangles in %s are not numbered correctly.\n",
1761 fname);
1762 printf(" (starting with triangle %d).\n", i);
1763 smallerr = 0;
1764 }
1765 stringptr = findfield(stringptr);
1766 if (*stringptr == '\0') {
1767 printf("Error: Triangle %d has no subdomain in %s.\n", i, fname);
1768 free(*partition);
1769 return 1;
1770 }
1771 (*partition)[i] = (int) strtol (stringptr, &stringptr, 0) - firstnumber;
1772 if (((*partition)[i] >= *parts) || ((*partition)[i] < 0)) {
1773 printf(" Error: Triangle %d of %s has an invalid subdomain.\n",
1774 i, fname);
1775 free(*partition);
1776 return 1;
1777 }
1778 }
1779 fclose(infile);
1780 *partcenter = (REAL *) malloc(((*parts + 1) << 1) * sizeof(REAL));
1781 if (*partcenter == (REAL *) NULL) {
1782 printf("Error: Out of memory.\n");
1783 free(*partition);
1784 return 1;
1785 }
1786 *partshift = (REAL *) malloc((*parts << 1) * sizeof(REAL));
1787 if (*partshift == (REAL *) NULL) {
1788 printf("Error: Out of memory.\n");
1789 free(*partition);
1790 free(*partcenter);
1791 return 1;
1792 }
1793 partsize = (int *) malloc((*parts + 1) * sizeof(int));
1794 if (partsize == (int *) NULL) {
1795 printf("Error: Out of memory.\n");
1796 free(*partition);
1797 free(*partcenter);
1798 free(*partshift);
1799 return 1;
1800 }
1801 index = 3;
1802 for (i = 0; i <= *parts; i++) {
1803 partsize[i] = 0;
1804 (*partcenter)[i << 1] = 0.0;
1805 (*partcenter)[(i << 1) + 1] = 0.0;
1806 }
1807 for (i = 1; i <= elems; i++) {
1808 partsize[(*partition)[i]] += 3;
1809 for (j = 0; j < 3; j++) {
1810 (*partcenter)[(*partition)[i] << 1] +=
1811 nodeptr[eleptr[index] * dim];
1812 (*partcenter)[((*partition)[i] << 1) + 1] +=
1813 nodeptr[eleptr[index++] * dim + 1];
1814 }
1815 }
1816 for (i = 0; i < *parts; i++) {
1817 (*partcenter)[i << 1] /= (REAL) partsize[i];
1818 (*partcenter)[(i << 1) + 1] /= (REAL) partsize[i];
1819 (*partcenter)[*parts << 1] += (*partcenter)[i << 1];
1820 (*partcenter)[(*parts << 1) + 1] += (*partcenter)[(i << 1) + 1];
1821 }
1822 (*partcenter)[*parts << 1] /= (REAL) *parts;
1823 (*partcenter)[(*parts << 1) + 1] /= (REAL) *parts;
1824 free(partsize);
1825 return 0;
1826 }
1827
load_adj(fname,subdomains,ptr)1828 int load_adj(fname, subdomains, ptr)
1829 char *fname;
1830 int *subdomains;
1831 int **ptr;
1832 {
1833 FILE *infile;
1834 char inputline[INPUTLINESIZE];
1835 char *stringptr;
1836 int i, j;
1837
1838 if (!quiet) {
1839 printf("Opening %s.\n", fname);
1840 }
1841 infile = fopen(fname, "r");
1842 if (infile == (FILE *) NULL) {
1843 printf(" Error: Cannot access file %s.\n", fname);
1844 return 1;
1845 }
1846 stringptr = readline(inputline, infile, fname);
1847 *subdomains = (int) strtol (stringptr, &stringptr, 0);
1848 if (*subdomains < 1) {
1849 printf(" Error: %s contains %d subdomains.\n", fname, *subdomains);
1850 return 1;
1851 }
1852 *ptr = (int *) malloc(*subdomains * *subdomains * sizeof(int));
1853 if (*ptr == (int *) NULL) {
1854 printf(" Out of memory.\n");
1855 return 1;
1856 }
1857 for (i = 0; i < *subdomains; i++) {
1858 for (j = 0; j < *subdomains; j++) {
1859 stringptr = readline(inputline, infile, fname);
1860 (*ptr)[i * *subdomains + j] = (int) strtol (stringptr, &stringptr, 0);
1861 }
1862 }
1863 return 0;
1864 }
1865
findpartshift(parts,explosion,partcenter,partshift)1866 void findpartshift(parts, explosion, partcenter, partshift)
1867 int parts;
1868 REAL explosion;
1869 REAL *partcenter;
1870 REAL *partshift;
1871 {
1872 int i;
1873
1874 for (i = 0; i < parts; i++) {
1875 partshift[i << 1] = explosion *
1876 (partcenter[i << 1] - partcenter[parts << 1]);
1877 partshift[(i << 1) + 1] = explosion *
1878 (partcenter[(i << 1) + 1] - partcenter[(parts << 1) + 1]);
1879 }
1880 }
1881
load_image(inc,image)1882 int load_image(inc, image)
1883 int inc;
1884 int image;
1885 {
1886 int error;
1887
1888 switch (image) {
1889 case NODE:
1890 error = load_node(nodefilename[inc], &firstnumber[inc], &nodes[inc],
1891 &node_dim[inc], &nodeptr[inc], &xlo[inc][NODE],
1892 &ylo[inc][NODE], &xhi[inc][NODE], &yhi[inc][NODE]);
1893 break;
1894 case POLY:
1895 error = load_poly(inc, polyfilename[inc], &firstnumber[inc],
1896 &polynodes[inc], &poly_dim[inc], &polyedges[inc],
1897 &polyholes[inc], &polynodeptr[inc], &polyedgeptr[inc],
1898 &polyholeptr[inc],
1899 &xlo[inc][POLY], &ylo[inc][POLY],
1900 &xhi[inc][POLY], &yhi[inc][POLY]);
1901 break;
1902 case ELE:
1903 error = load_ele(elefilename[inc], firstnumber[inc], nodes[inc],
1904 &elems[inc], &ele_corners[inc], &eleptr[inc]);
1905 xlo[inc][ELE] = xlo[inc][NODE];
1906 ylo[inc][ELE] = ylo[inc][NODE];
1907 xhi[inc][ELE] = xhi[inc][NODE];
1908 yhi[inc][ELE] = yhi[inc][NODE];
1909 break;
1910 case EDGE:
1911 error = load_edge(edgefilename[inc], firstnumber[inc], nodes[inc],
1912 &edges[inc], &edgeptr[inc], &normptr[inc]);
1913 xlo[inc][EDGE] = xlo[inc][NODE];
1914 ylo[inc][EDGE] = ylo[inc][NODE];
1915 xhi[inc][EDGE] = xhi[inc][NODE];
1916 yhi[inc][EDGE] = yhi[inc][NODE];
1917 break;
1918 case PART:
1919 error = load_part(partfilename[inc], node_dim[inc], firstnumber[inc],
1920 elems[inc], nodeptr[inc], eleptr[inc],
1921 &subdomains[inc], &partpart[inc], &partcenter[inc],
1922 &partshift[inc]);
1923 if (!error) {
1924 findpartshift(subdomains[inc], explosion, partcenter[inc],
1925 partshift[inc]);
1926 }
1927 xlo[inc][PART] = xlo[inc][NODE];
1928 ylo[inc][PART] = ylo[inc][NODE];
1929 xhi[inc][PART] = xhi[inc][NODE];
1930 yhi[inc][PART] = yhi[inc][NODE];
1931 break;
1932 case ADJ:
1933 error = load_adj(adjfilename[inc], &adjsubdomains[inc], &adjptr[inc]);
1934 xlo[inc][ADJ] = xlo[inc][NODE];
1935 ylo[inc][ADJ] = ylo[inc][NODE];
1936 xhi[inc][ADJ] = xhi[inc][NODE];
1937 yhi[inc][ADJ] = yhi[inc][NODE];
1938 break;
1939 case VORO:
1940 error = load_node(vnodefilename[inc], &firstnumber[inc], &vnodes[inc],
1941 &vnode_dim[inc], &vnodeptr[inc], &xlo[inc][VORO],
1942 &ylo[inc][VORO], &xhi[inc][VORO], &yhi[inc][VORO]);
1943 if (!error) {
1944 error = load_edge(vedgefilename[inc], firstnumber[inc], vnodes[inc],
1945 &vedges[inc], &vedgeptr[inc], &vnormptr[inc]);
1946 }
1947 break;
1948 default:
1949 error = 1;
1950 }
1951 if (!error) {
1952 loaded[inc][image] = 1;
1953 }
1954 return error;
1955 }
1956
choose_image(inc,image)1957 void choose_image(inc, image)
1958 int inc;
1959 int image;
1960 {
1961 if (!loaded[inc][image]) {
1962 if ((image == ELE) || (image == EDGE) || (image == PART)
1963 || (image == ADJ)) {
1964 if (!loaded[inc][NODE]) {
1965 if (load_image(inc, NODE)) {
1966 return;
1967 }
1968 }
1969 }
1970 if ((image == PART) || (image == ADJ)) {
1971 if (!loaded[inc][ELE]) {
1972 if (load_image(inc, ELE)) {
1973 return;
1974 }
1975 }
1976 }
1977 if (image == ADJ) {
1978 if (!loaded[inc][PART]) {
1979 if (load_image(inc, PART)) {
1980 return;
1981 }
1982 }
1983 }
1984 if (load_image(inc, image)) {
1985 return;
1986 }
1987 }
1988 current_inc = inc;
1989 current_image = image;
1990 }
1991
make_button(name,x,y,width)1992 Window make_button(name, x, y, width)
1993 char *name;
1994 int x;
1995 int y;
1996 int width;
1997 {
1998 XSetWindowAttributes attr;
1999 XSizeHints hints;
2000 Window button;
2001
2002 attr.background_pixel = black;
2003 attr.border_pixel = white;
2004 attr.backing_store = NotUseful;
2005 attr.event_mask = ExposureMask | ButtonReleaseMask | ButtonPressMask;
2006 attr.bit_gravity = SouthWestGravity;
2007 attr.win_gravity = SouthWestGravity;
2008 attr.save_under = False;
2009 button = XCreateWindow(display, mainwindow, x, y, width, BUTTONHEIGHT - 4,
2010 2, 0, InputOutput, CopyFromParent,
2011 CWBackPixel | CWBorderPixel | CWEventMask |
2012 CWBitGravity | CWWinGravity | CWBackingStore |
2013 CWSaveUnder, &attr);
2014 hints.width = width;
2015 hints.height = BUTTONHEIGHT - 4;
2016 hints.min_width = 0;
2017 hints.min_height = BUTTONHEIGHT - 4;
2018 hints.max_width = width;
2019 hints.max_height = BUTTONHEIGHT - 4;
2020 hints.width_inc = 1;
2021 hints.height_inc = 1;
2022 hints.flags = PMinSize | PMaxSize | PSize | PResizeInc;
2023 XSetStandardProperties(display, button, name, "showme", None, (char **) NULL,
2024 0, &hints);
2025 return button;
2026 }
2027
make_buttons(y)2028 void make_buttons(y)
2029 int y;
2030 {
2031 int i;
2032
2033 for (i = 1; i >= 0; i--) {
2034 nodewin[i] = make_button("node", 0, y + (1 - i) * BUTTONHEIGHT, 42);
2035 XMapWindow(display, nodewin[i]);
2036 polywin[i] = make_button("poly", 44, y + (1 - i) * BUTTONHEIGHT, 42);
2037 XMapWindow(display, polywin[i]);
2038 elewin[i] = make_button("ele", 88, y + (1 - i) * BUTTONHEIGHT, 33);
2039 XMapWindow(display, elewin[i]);
2040 edgewin[i] = make_button("edge", 123, y + (1 - i) * BUTTONHEIGHT, 42);
2041 XMapWindow(display, edgewin[i]);
2042 partwin[i] = make_button("part", 167, y + (1 - i) * BUTTONHEIGHT, 42);
2043 XMapWindow(display, partwin[i]);
2044 adjwin[i] = make_button("adj", 211, y + (1 - i) * BUTTONHEIGHT, 33);
2045 XMapWindow(display, adjwin[i]);
2046 voronoiwin[i] = make_button("voro", 246, y + (1 - i) * BUTTONHEIGHT, 42);
2047 XMapWindow(display, voronoiwin[i]);
2048 }
2049 versionpluswin = make_button(" +", 290, y, 52);
2050 XMapWindow(display, versionpluswin);
2051 versionminuswin = make_button(" -", 290, y + BUTTONHEIGHT, 52);
2052 XMapWindow(display, versionminuswin);
2053
2054 quitwin = make_button("Quit", 0, y + 2 * BUTTONHEIGHT, 42);
2055 XMapWindow(display, quitwin);
2056 leftwin = make_button("<", 44, y + 2 * BUTTONHEIGHT, 14);
2057 XMapWindow(display, leftwin);
2058 rightwin = make_button(">", 60, y + 2 * BUTTONHEIGHT, 14);
2059 XMapWindow(display, rightwin);
2060 upwin = make_button("^", 76, y + 2 * BUTTONHEIGHT, 14);
2061 XMapWindow(display, upwin);
2062 downwin = make_button("v", 92, y + 2 * BUTTONHEIGHT, 14);
2063 XMapWindow(display, downwin);
2064 resetwin = make_button("Reset", 108, y + 2 * BUTTONHEIGHT, 52);
2065 XMapWindow(display, resetwin);
2066 widthpluswin = make_button("Width+", 162, y + 2 * BUTTONHEIGHT, 61);
2067 XMapWindow(display, widthpluswin);
2068 widthminuswin = make_button("-", 225, y + 2 * BUTTONHEIGHT, 14);
2069 XMapWindow(display, widthminuswin);
2070 expwin = make_button("Exp", 241, y + 2 * BUTTONHEIGHT, 33);
2071 XMapWindow(display, expwin);
2072 exppluswin = make_button("+", 276, y + 2 * BUTTONHEIGHT, 14);
2073 XMapWindow(display, exppluswin);
2074 expminuswin = make_button("-", 292, y + 2 * BUTTONHEIGHT, 14);
2075 XMapWindow(display, expminuswin);
2076 fillwin = make_button("Fill", 308, y + 2 * BUTTONHEIGHT, 41);
2077 XMapWindow(display, fillwin);
2078 pswin = make_button("PS", 351, y + 2 * BUTTONHEIGHT, 24);
2079 XMapWindow(display, pswin);
2080 epswin = make_button("EPS", 377, y + 2 * BUTTONHEIGHT, 33);
2081 XMapWindow(display, epswin);
2082 }
2083
fill_button(button)2084 void fill_button(button)
2085 Window button;
2086 {
2087 int x, y;
2088 unsigned int w, h, d, b;
2089 Window rootw;
2090
2091 XGetGeometry(display, button, &rootw, &x, &y, &w, &h, &d, &b);
2092 XFillRectangle(display, button, fontgc, 0, 0, w, h);
2093 }
2094
draw_buttons()2095 void draw_buttons()
2096 {
2097 char numberstring[32];
2098 char buttonstring[6];
2099 int i;
2100
2101 for (i = 1; i >= 0; i--) {
2102 if ((current_image == NODE) && (current_inc == i)) {
2103 fill_button(nodewin[i]);
2104 XDrawString(display, nodewin[i], blackfontgc, 2, 13, "node", 4);
2105 } else {
2106 XClearWindow(display, nodewin[i]);
2107 XDrawString(display, nodewin[i], fontgc, 2, 13, "node", 4);
2108 }
2109 if ((current_image == POLY) && (current_inc == i)) {
2110 fill_button(polywin[i]);
2111 XDrawString(display, polywin[i], blackfontgc, 2, 13, "poly", 4);
2112 } else {
2113 XClearWindow(display, polywin[i]);
2114 XDrawString(display, polywin[i], fontgc, 2, 13, "poly", 4);
2115 }
2116 if ((current_image == ELE) && (current_inc == i)) {
2117 fill_button(elewin[i]);
2118 XDrawString(display, elewin[i], blackfontgc, 2, 13, "ele", 3);
2119 } else {
2120 XClearWindow(display, elewin[i]);
2121 XDrawString(display, elewin[i], fontgc, 2, 13, "ele", 3);
2122 }
2123 if ((current_image == EDGE) && (current_inc == i)) {
2124 fill_button(edgewin[i]);
2125 XDrawString(display, edgewin[i], blackfontgc, 2, 13, "edge", 4);
2126 } else {
2127 XClearWindow(display, edgewin[i]);
2128 XDrawString(display, edgewin[i], fontgc, 2, 13, "edge", 4);
2129 }
2130 if ((current_image == PART) && (current_inc == i)) {
2131 fill_button(partwin[i]);
2132 XDrawString(display, partwin[i], blackfontgc, 2, 13, "part", 4);
2133 } else {
2134 XClearWindow(display, partwin[i]);
2135 XDrawString(display, partwin[i], fontgc, 2, 13, "part", 4);
2136 }
2137 if ((current_image == ADJ) && (current_inc == i)) {
2138 fill_button(adjwin[i]);
2139 XDrawString(display, adjwin[i], blackfontgc, 2, 13, "adj", 3);
2140 } else {
2141 XClearWindow(display, adjwin[i]);
2142 XDrawString(display, adjwin[i], fontgc, 2, 13, "adj", 3);
2143 }
2144 if ((current_image == VORO) && (current_inc == i)) {
2145 fill_button(voronoiwin[i]);
2146 XDrawString(display, voronoiwin[i], blackfontgc, 2, 13, "voro", 4);
2147 } else {
2148 XClearWindow(display, voronoiwin[i]);
2149 XDrawString(display, voronoiwin[i], fontgc, 2, 13, "voro", 4);
2150 }
2151 }
2152
2153 XClearWindow(display, versionpluswin);
2154 sprintf(numberstring, "%d", loweriteration + 1);
2155 sprintf(buttonstring, "%-4.4s+", numberstring);
2156 XDrawString(display, versionpluswin, fontgc, 2, 13, buttonstring, 5);
2157 XClearWindow(display, versionminuswin);
2158 sprintf(numberstring, "%d", loweriteration);
2159 if (loweriteration == 0) {
2160 sprintf(buttonstring, "%-4.4s", numberstring);
2161 } else {
2162 sprintf(buttonstring, "%-4.4s-", numberstring);
2163 }
2164 XDrawString(display, versionminuswin, fontgc, 2, 13, buttonstring, 5);
2165
2166 XClearWindow(display, quitwin);
2167 XDrawString(display, quitwin, fontgc, 2, 13, "Quit", 4);
2168 XClearWindow(display, leftwin);
2169 XDrawString(display, leftwin, fontgc, 2, 13, "<", 1);
2170 XClearWindow(display, rightwin);
2171 XDrawString(display, rightwin, fontgc, 2, 13, ">", 1);
2172 XClearWindow(display, upwin);
2173 XDrawString(display, upwin, fontgc, 2, 13, "^", 1);
2174 XClearWindow(display, downwin);
2175 XDrawString(display, downwin, fontgc, 2, 13, "v", 1);
2176 XClearWindow(display, resetwin);
2177 XDrawString(display, resetwin, fontgc, 2, 13, "Reset", 6);
2178 XClearWindow(display, widthpluswin);
2179 if (line_width < 100) {
2180 XDrawString(display, widthpluswin, fontgc, 2, 13, "Width+", 6);
2181 } else {
2182 XDrawString(display, widthpluswin, fontgc, 2, 13, "Width ", 6);
2183 }
2184 XClearWindow(display, widthminuswin);
2185 if (line_width > 1) {
2186 XDrawString(display, widthminuswin, fontgc, 2, 13, "-", 1);
2187 }
2188 XClearWindow(display, expwin);
2189 XClearWindow(display, exppluswin);
2190 XClearWindow(display, expminuswin);
2191 XClearWindow(display, fillwin);
2192 if (current_image == PART) {
2193 if (explode) {
2194 fill_button(expwin);
2195 XDrawString(display, expwin, blackfontgc, 2, 13, "Exp", 3);
2196 } else {
2197 XDrawString(display, expwin, fontgc, 2, 13, "Exp", 3);
2198 }
2199 XDrawString(display, exppluswin, fontgc, 2, 13, "+", 1);
2200 XDrawString(display, expminuswin, fontgc, 2, 13, "-", 1);
2201 if (fillelem) {
2202 fill_button(fillwin);
2203 XDrawString(display, fillwin, blackfontgc, 2, 13, "Fill", 4);
2204 } else {
2205 XDrawString(display, fillwin, fontgc, 2, 13, "Fill", 4);
2206 }
2207 }
2208 XClearWindow(display, pswin);
2209 XDrawString(display, pswin, fontgc, 2, 13, "PS", 2);
2210 XClearWindow(display, epswin);
2211 XDrawString(display, epswin, fontgc, 2, 13, "EPS", 3);
2212 }
2213
showme_window(argc,argv)2214 void showme_window(argc, argv)
2215 int argc;
2216 char **argv;
2217 {
2218 XSetWindowAttributes attr;
2219 XSizeHints hints;
2220 XGCValues fontvalues, linevalues;
2221 XColor alloc_color, exact_color;
2222 int i;
2223
2224 display = XOpenDisplay((char *) NULL);
2225 if (!display) {
2226 printf("Error: Cannot open display.\n");
2227 exit(1);
2228 }
2229 screen = DefaultScreen(display);
2230 rootwindow = DefaultRootWindow(display);
2231 black = BlackPixel(display, screen);
2232 white = WhitePixel(display, screen);
2233 windowdepth = DefaultDepth(display, screen);
2234 rootmap = DefaultColormap(display, screen);
2235 width = STARTWIDTH;
2236 height = STARTHEIGHT;
2237 attr.background_pixel = black;
2238 attr.border_pixel = white;
2239 attr.backing_store = NotUseful;
2240 attr.event_mask = ExposureMask | ButtonReleaseMask | ButtonPressMask |
2241 StructureNotifyMask;
2242 attr.bit_gravity = NorthWestGravity;
2243 attr.win_gravity = NorthWestGravity;
2244 attr.save_under = False;
2245 mainwindow = XCreateWindow(display, rootwindow, 0, 0, width,
2246 height + PANELHEIGHT, 3, 0,
2247 InputOutput, CopyFromParent,
2248 CWBackPixel | CWBorderPixel | CWEventMask |
2249 CWBitGravity | CWWinGravity | CWBackingStore |
2250 CWSaveUnder, &attr);
2251 hints.width = width;
2252 hints.height = height + PANELHEIGHT;
2253 hints.min_width = MINWIDTH;
2254 hints.min_height = MINHEIGHT + PANELHEIGHT;
2255 hints.width_inc = 1;
2256 hints.height_inc = 1;
2257 hints.flags = PMinSize | PSize | PResizeInc;
2258 XSetStandardProperties(display, mainwindow, "Show Me", "showme", None,
2259 argv, argc, &hints);
2260
2261 static const unsigned char temp_show_me_achimedes_local[18] = {'s','h','o','w','m','e','\0','A','r','c','h','i','m','e','d','e','s','\0'};
2262 XChangeProperty(display, mainwindow, XA_WM_CLASS, XA_STRING, 8,
2263 PropModeReplace, temp_show_me_achimedes_local, 18U);
2264 XClearWindow(display, mainwindow);
2265 XMapWindow(display, mainwindow);
2266 if ((windowdepth > 1) &&
2267 XAllocNamedColor(display, rootmap, "yellow", &alloc_color,
2268 &exact_color)) {
2269 color = 1;
2270 explode = bw_ps;
2271 fontvalues.foreground = alloc_color.pixel;
2272 linevalues.foreground = alloc_color.pixel;
2273 showme_foreground = alloc_color.pixel;
2274 for (i = 0; i < 64; i++) {
2275 if (XAllocNamedColor(display, rootmap, colorname[i], &alloc_color,
2276 &rgb[i])) {
2277 colors[i] = alloc_color.pixel;
2278 } else {
2279 colors[i] = white;
2280 rgb[i].red = alloc_color.red;
2281 rgb[i].green = alloc_color.green;
2282 rgb[i].blue = alloc_color.blue;
2283 if (!quiet) {
2284 printf("Warning: I could not allocate %s.\n", colorname[i]);
2285 }
2286 }
2287 }
2288 } else {
2289 color = 0;
2290 fillelem = 0;
2291 explode = 1;
2292 fontvalues.foreground = white;
2293 linevalues.foreground = white;
2294 showme_foreground = white;
2295 }
2296 font = XLoadQueryFont(display, "9x15");
2297 fontvalues.background = black;
2298 fontvalues.font = font->fid;
2299 fontvalues.fill_style = FillSolid;
2300 fontvalues.line_width = 2;
2301 fontgc = XCreateGC(display, rootwindow, GCForeground | GCBackground |
2302 GCFont | GCLineWidth | GCFillStyle, &fontvalues);
2303 fontvalues.foreground = black;
2304 blackfontgc = XCreateGC(display, rootwindow, GCForeground | GCBackground |
2305 GCFont | GCLineWidth | GCFillStyle, &fontvalues);
2306 linevalues.background = black;
2307 linevalues.line_width = line_width;
2308 linevalues.cap_style = CapRound;
2309 linevalues.join_style = JoinRound;
2310 linevalues.fill_style = FillSolid;
2311 linegc = XCreateGC(display, rootwindow, GCForeground | GCBackground |
2312 GCLineWidth | GCCapStyle | GCJoinStyle | GCFillStyle,
2313 &linevalues);
2314 trianglegc = XCreateGC(display, rootwindow, GCForeground | GCBackground |
2315 GCLineWidth | GCCapStyle | GCJoinStyle | GCFillStyle,
2316 &linevalues);
2317 make_buttons(height);
2318 XFlush(display);
2319 }
2320
draw_node(nodes,dim,ptr,xscale,yscale,xoffset,yoffset)2321 void draw_node(nodes, dim, ptr, xscale, yscale, xoffset, yoffset)
2322 int nodes;
2323 int dim;
2324 REAL *ptr;
2325 REAL xscale;
2326 REAL yscale;
2327 REAL xoffset;
2328 REAL yoffset;
2329 {
2330 int i;
2331 int index;
2332
2333 index = dim;
2334 for (i = 1; i <= nodes; i++) {
2335 XFillRectangle(display, mainwindow, linegc,
2336 (int) (ptr[index] * xscale + xoffset) - (line_width >> 1),
2337 (int) (ptr[index + 1] * yscale + yoffset) -
2338 (line_width >> 1), line_width, line_width);
2339 index += dim;
2340 }
2341 }
2342
draw_poly(nodes,dim,edges,holes,nodeptr,edgeptr,holeptr,xscale,yscale,xoffset,yoffset)2343 void draw_poly(nodes, dim, edges, holes, nodeptr, edgeptr, holeptr,
2344 xscale, yscale, xoffset, yoffset)
2345 int nodes;
2346 int dim;
2347 int edges;
2348 int holes;
2349 REAL *nodeptr;
2350 int *edgeptr;
2351 REAL *holeptr;
2352 REAL xscale;
2353 REAL yscale;
2354 REAL xoffset;
2355 REAL yoffset;
2356 {
2357 int i;
2358 int index;
2359 REAL *point1, *point2;
2360 int x1, y1, x2, y2;
2361
2362 index = dim;
2363 for (i = 1; i <= nodes; i++) {
2364 XFillRectangle(display, mainwindow, linegc,
2365 (int) (nodeptr[index] * xscale + xoffset) -
2366 (line_width >> 1),
2367 (int) (nodeptr[index + 1] * yscale + yoffset) -
2368 (line_width >> 1), line_width, line_width);
2369 index += dim;
2370 }
2371 index = 2;
2372 for (i = 1; i <= edges; i++) {
2373 point1 = &nodeptr[edgeptr[index++] * dim];
2374 point2 = &nodeptr[edgeptr[index++] * dim];
2375 XDrawLine(display, mainwindow, linegc,
2376 (int) (point1[0] * xscale + xoffset),
2377 (int) (point1[1] * yscale + yoffset),
2378 (int) (point2[0] * xscale + xoffset),
2379 (int) (point2[1] * yscale + yoffset));
2380 }
2381 index = dim;
2382 if (color) {
2383 XSetForeground(display, linegc, colors[0]);
2384 }
2385 for (i = 1; i <= holes; i++) {
2386 x1 = (int) (holeptr[index] * xscale + xoffset) - 3;
2387 y1 = (int) (holeptr[index + 1] * yscale + yoffset) - 3;
2388 x2 = x1 + 6;
2389 y2 = y1 + 6;
2390 XDrawLine(display, mainwindow, linegc, x1, y1, x2, y2);
2391 XDrawLine(display, mainwindow, linegc, x1, y2, x2, y1);
2392 index += dim;
2393 }
2394 XSetForeground(display, linegc, showme_foreground);
2395 }
2396
draw_ele(inc,elems,corners,ptr,partition,shift,xscale,yscale,xoffset,yoffset)2397 void draw_ele(inc, elems, corners, ptr, partition, shift,
2398 xscale, yscale, xoffset, yoffset)
2399 int inc;
2400 int elems;
2401 int corners; /* unused */
2402 int *ptr;
2403 int *partition;
2404 REAL *shift;
2405 REAL xscale;
2406 REAL yscale;
2407 REAL xoffset;
2408 REAL yoffset;
2409 {
2410 int i, j;
2411 int index;
2412 REAL shiftx = 0.0, shifty = 0.0;
2413 REAL *prevpoint, *nowpoint;
2414 XPoint *vertices = (XPoint *) NULL;
2415
2416 if (color && fillelem && (partition != (int *) NULL)) {
2417 vertices = (XPoint *) malloc(3 * sizeof(XPoint));
2418 if (vertices == (XPoint *) NULL) {
2419 printf("Error: Out of memory.\n");
2420 exit(1);
2421 }
2422 }
2423 index = 3;
2424 for (i = 1; i <= elems; i++) {
2425 if ((partition != (int *) NULL) && explode) {
2426 shiftx = shift[partition[i] << 1];
2427 shifty = shift[(partition[i] << 1) + 1];
2428 }
2429 if (color && (partition != (int *) NULL)) {
2430 if (fillelem) {
2431 XSetForeground(display, trianglegc, colors[partition[i] & 63]);
2432 } else {
2433 XSetForeground(display, linegc, colors[partition[i] & 63]);
2434 }
2435 }
2436 if (color && fillelem && (partition != (int *) NULL)) {
2437 if ((partition != (int *) NULL) && explode) {
2438 for (j = 0; j < 3; j++) {
2439 nowpoint = &nodeptr[inc][ptr[index + j] * node_dim[inc]];
2440 vertices[j].x = (nowpoint[0] + shiftx) * xscale + xoffset;
2441 vertices[j].y = (nowpoint[1] + shifty) * yscale + yoffset;
2442 }
2443 } else {
2444 for (j = 0; j < 3; j++) {
2445 nowpoint = &nodeptr[inc][ptr[index + j] * node_dim[inc]];
2446 vertices[j].x = nowpoint[0] * xscale + xoffset;
2447 vertices[j].y = nowpoint[1] * yscale + yoffset;
2448 }
2449 }
2450 XFillPolygon(display, mainwindow, trianglegc, vertices, 3,
2451 Convex, CoordModeOrigin);
2452 }
2453 prevpoint = &nodeptr[inc][ptr[index + 2] * node_dim[inc]];
2454 if ((partition != (int *) NULL) && explode) {
2455 for (j = 0; j < 3; j++) {
2456 nowpoint = &nodeptr[inc][ptr[index++] * node_dim[inc]];
2457 XDrawLine(display, mainwindow, linegc,
2458 (int) ((prevpoint[0] + shiftx) * xscale + xoffset),
2459 (int) ((prevpoint[1] + shifty) * yscale + yoffset),
2460 (int) ((nowpoint[0] + shiftx) * xscale + xoffset),
2461 (int) ((nowpoint[1] + shifty) * yscale + yoffset));
2462 prevpoint = nowpoint;
2463 }
2464 } else {
2465 for (j = 0; j < 3; j++) {
2466 nowpoint = &nodeptr[inc][ptr[index++] * node_dim[inc]];
2467 XDrawLine(display, mainwindow, linegc,
2468 (int) (prevpoint[0] * xscale + xoffset),
2469 (int) (prevpoint[1] * yscale + yoffset),
2470 (int) (nowpoint[0] * xscale + xoffset),
2471 (int) (nowpoint[1] * yscale + yoffset));
2472 prevpoint = nowpoint;
2473 }
2474 }
2475 }
2476 if (color && fillelem && (partition != (int *) NULL)) {
2477 free(vertices);
2478 }
2479 XSetForeground(display, linegc, showme_foreground);
2480 }
2481
draw_edge(nodes,dim,edges,nodeptr,edgeptr,normptr,xscale,yscale,xoffset,yoffset)2482 void draw_edge(nodes, dim, edges, nodeptr, edgeptr, normptr,
2483 xscale, yscale, xoffset, yoffset)
2484 int nodes; /* unused */
2485 int dim;
2486 int edges;
2487 REAL *nodeptr;
2488 int *edgeptr;
2489 REAL *normptr;
2490 REAL xscale;
2491 REAL yscale;
2492 REAL xoffset;
2493 REAL yoffset;
2494 {
2495 int i;
2496 int index;
2497 REAL *point1, *point2;
2498 REAL normx, normy;
2499 REAL normmult, normmultx, normmulty;
2500 REAL windowxmin, windowymin, windowxmax, windowymax;
2501
2502 index = 2;
2503 for (i = 1; i <= edges; i++) {
2504 point1 = &nodeptr[edgeptr[index++] * dim];
2505 if (edgeptr[index] == -1) {
2506 normx = normptr[index - 1];
2507 normy = normptr[index++];
2508 normmultx = 0.0;
2509 if (normx > 0) {
2510 windowxmax = (width - 1 - xoffset) / xscale;
2511 normmultx = (windowxmax - point1[0]) / normx;
2512 } else if (normx < 0) {
2513 windowxmin = -xoffset / xscale;
2514 normmultx = (windowxmin - point1[0]) / normx;
2515 }
2516 normmulty = 0.0;
2517 if (normy > 0) {
2518 windowymax = -yoffset / yscale;
2519 normmulty = (windowymax - point1[1]) / normy;
2520 } else if (normy < 0) {
2521 windowymin = (height - 1 - yoffset) / yscale;
2522 normmulty = (windowymin - point1[1]) / normy;
2523 }
2524 if (normmultx == 0.0) {
2525 normmult = normmulty;
2526 } else if (normmulty == 0.0) {
2527 normmult = normmultx;
2528 } else if (normmultx < normmulty) {
2529 normmult = normmultx;
2530 } else {
2531 normmult = normmulty;
2532 }
2533 if (normmult > 0.0) {
2534 XDrawLine(display, mainwindow, linegc,
2535 (int) (point1[0] * xscale + xoffset),
2536 (int) (point1[1] * yscale + yoffset),
2537 (int) ((point1[0] + normmult * normx) * xscale + xoffset),
2538 (int) ((point1[1] + normmult * normy) * yscale + yoffset));
2539 }
2540 } else {
2541 point2 = &nodeptr[edgeptr[index++] * dim];
2542 XDrawLine(display, mainwindow, linegc,
2543 (int) (point1[0] * xscale + xoffset),
2544 (int) (point1[1] * yscale + yoffset),
2545 (int) (point2[0] * xscale + xoffset),
2546 (int) (point2[1] * yscale + yoffset));
2547 }
2548 }
2549 }
2550
draw_adj(dim,subdomains,ptr,center,xscale,yscale,xoffset,yoffset)2551 void draw_adj(dim, subdomains, ptr, center, xscale, yscale,
2552 xoffset, yoffset)
2553 int dim;
2554 int subdomains;
2555 int *ptr;
2556 REAL *center;
2557 REAL xscale;
2558 REAL yscale;
2559 REAL xoffset;
2560 REAL yoffset;
2561 {
2562 int i, j;
2563 REAL *point1, *point2;
2564
2565 for (i = 0; i < subdomains; i++) {
2566 for (j = i + 1; j < subdomains; j++) {
2567 if (ptr[i * subdomains + j]) {
2568 point1 = ¢er[i * dim];
2569 point2 = ¢er[j * dim];
2570 XDrawLine(display, mainwindow, linegc,
2571 (int) (point1[0] * xscale + xoffset),
2572 (int) (point1[1] * yscale + yoffset),
2573 (int) (point2[0] * xscale + xoffset),
2574 (int) (point2[1] * yscale + yoffset));
2575 }
2576 }
2577 }
2578 for (i = 0; i < subdomains; i++) {
2579 point1 = ¢er[i * dim];
2580 if (color) {
2581 XSetForeground(display, linegc, colors[i & 63]);
2582 }
2583 XFillArc(display, mainwindow, linegc,
2584 (int) (point1[0] * xscale + xoffset) - 5 - (line_width >> 1),
2585 (int) (point1[1] * yscale + yoffset) - 5 - (line_width >> 1),
2586 line_width + 10, line_width + 10, 0, 23040);
2587 }
2588 XSetForeground(display, linegc, showme_foreground);
2589 }
2590
draw(inc,image,xmin,ymin,xmax,ymax)2591 void draw(inc, image, xmin, ymin, xmax, ymax)
2592 int inc;
2593 int image;
2594 REAL xmin;
2595 REAL ymin;
2596 REAL xmax;
2597 REAL ymax;
2598 {
2599 draw_buttons();
2600 XClearWindow(display, mainwindow);
2601 if (image == NOTHING) {
2602 return;
2603 }
2604 if (!loaded[inc][image]) {
2605 return;
2606 }
2607 if ((image == PART) && explode) {
2608 xmin += (xmin - partcenter[inc][subdomains[inc] << 1]) * explosion;
2609 xmax += (xmax - partcenter[inc][subdomains[inc] << 1]) * explosion;
2610 ymin += (ymin - partcenter[inc][(subdomains[inc] << 1) + 1]) * explosion;
2611 ymax += (ymax - partcenter[inc][(subdomains[inc] << 1) + 1]) * explosion;
2612 }
2613 xscale = (REAL) (width - line_width - 4) / (xmax - xmin);
2614 yscale = (REAL) (height - line_width - 4) / (ymax - ymin);
2615 if (xscale > yscale) {
2616 xscale = yscale;
2617 } else {
2618 yscale = xscale;
2619 }
2620 xoffset = 0.5 * ((REAL) width - xscale * (xmax - xmin)) -
2621 xscale * xmin;
2622 yoffset = (REAL) height - 0.5 * ((REAL) height - yscale * (ymax - ymin)) +
2623 yscale * ymin;
2624 yscale = - yscale;
2625 switch(image) {
2626 case NODE:
2627 draw_node(nodes[inc], node_dim[inc], nodeptr[inc],
2628 xscale, yscale, xoffset, yoffset);
2629 break;
2630 case POLY:
2631 if (polynodes[inc] > 0) {
2632 draw_poly(polynodes[inc], poly_dim[inc], polyedges[inc],
2633 polyholes[inc], polynodeptr[inc], polyedgeptr[inc],
2634 polyholeptr[inc],
2635 xscale, yscale, xoffset, yoffset);
2636 } else {
2637 draw_poly(nodes[inc], node_dim[inc], polyedges[inc],
2638 polyholes[inc], nodeptr[inc], polyedgeptr[inc],
2639 polyholeptr[inc],
2640 xscale, yscale, xoffset, yoffset);
2641 }
2642 break;
2643 case ELE:
2644 draw_ele(inc, elems[inc], ele_corners[inc], eleptr[inc],
2645 (int *) NULL, (REAL *) NULL,
2646 xscale, yscale, xoffset, yoffset);
2647 break;
2648 case EDGE:
2649 draw_edge(nodes[inc], node_dim[inc], edges[inc],
2650 nodeptr[inc], edgeptr[inc], normptr[inc],
2651 xscale, yscale, xoffset, yoffset);
2652 break;
2653 case PART:
2654 draw_ele(inc, elems[inc], ele_corners[inc], eleptr[inc],
2655 partpart[inc], partshift[inc],
2656 xscale, yscale, xoffset, yoffset);
2657 break;
2658 case ADJ:
2659 draw_adj(node_dim[inc], adjsubdomains[inc], adjptr[inc], partcenter[inc],
2660 xscale, yscale, xoffset, yoffset);
2661 break;
2662 case VORO:
2663 if (loaded[inc][NODE]) {
2664 draw_node(nodes[inc], node_dim[inc], nodeptr[inc],
2665 xscale, yscale, xoffset, yoffset);
2666 }
2667 draw_edge(vnodes[inc], vnode_dim[inc], vedges[inc],
2668 vnodeptr[inc], vedgeptr[inc], vnormptr[inc],
2669 xscale, yscale, xoffset, yoffset);
2670 break;
2671 default:
2672 break;
2673 }
2674 }
2675
addps(instring,outstring,eps)2676 void addps(instring, outstring, eps)
2677 char *instring;
2678 char *outstring;
2679 int eps;
2680 {
2681 strcpy(outstring, instring);
2682 if (eps) {
2683 strcat(outstring, ".eps");
2684 } else {
2685 strcat(outstring, ".ps");
2686 }
2687 }
2688
print_head(fname,file,llcornerx,llcornery,eps)2689 int print_head(fname, file, llcornerx, llcornery, eps)
2690 char *fname;
2691 FILE **file;
2692 int llcornerx;
2693 int llcornery;
2694 int eps;
2695 {
2696 if (!quiet) {
2697 printf("Writing %s\n", fname);
2698 }
2699 *file = fopen(fname, "w");
2700 if (*file == (FILE *) NULL) {
2701 printf(" Error: Could not open %s\n", fname);
2702 return 1;
2703 }
2704 if (eps) {
2705 fprintf(*file, "%%!PS-Adobe-2.0 EPSF-2.0\n");
2706 } else {
2707 fprintf(*file, "%%!PS-Adobe-2.0\n");
2708 }
2709 fprintf(*file, "%%%%BoundingBox: %d %d %d %d\n", llcornerx, llcornery,
2710 612 - llcornerx, 792 - llcornery);
2711 fprintf(*file, "%%%%Creator: Show Me\n");
2712 fprintf(*file, "%%%%EndComments\n\n");
2713 fprintf(*file, "1 setlinecap\n");
2714 fprintf(*file, "1 setlinejoin\n");
2715 fprintf(*file, "%d setlinewidth\n", line_width);
2716 fprintf(*file, "%d %d moveto\n", llcornerx, llcornery);
2717 fprintf(*file, "%d %d lineto\n", 612 - llcornerx, llcornery);
2718 fprintf(*file, "%d %d lineto\n", 612 - llcornerx, 792 - llcornery);
2719 fprintf(*file, "%d %d lineto\n", llcornerx, 792 - llcornery);
2720 fprintf(*file, "closepath\nclip\nnewpath\n");
2721 return 0;
2722 }
2723
print_node(nodefile,nodes,dim,ptr,xscale,yscale,xoffset,yoffset)2724 void print_node(nodefile, nodes, dim, ptr, xscale, yscale,
2725 xoffset, yoffset)
2726 FILE *nodefile;
2727 int nodes;
2728 int dim;
2729 REAL *ptr;
2730 REAL xscale;
2731 REAL yscale;
2732 REAL xoffset;
2733 REAL yoffset;
2734 {
2735 int i;
2736 int index;
2737
2738 index = dim;
2739 for (i = 1; i <= nodes; i++) {
2740 fprintf(nodefile, "%d %d %d 0 360 arc\nfill\n",
2741 (int) (ptr[index] * xscale + xoffset),
2742 (int) (ptr[index + 1] * yscale + yoffset),
2743 1 + (line_width >> 1));
2744 index += dim;
2745 }
2746 }
2747
print_poly(polyfile,nodes,dim,edges,holes,nodeptr,edgeptr,holeptr,xscale,yscale,xoffset,yoffset)2748 void print_poly(polyfile, nodes, dim, edges, holes, nodeptr, edgeptr, holeptr,
2749 xscale, yscale, xoffset, yoffset)
2750 FILE *polyfile;
2751 int nodes;
2752 int dim;
2753 int edges;
2754 int holes; /* unused */
2755 REAL *nodeptr;
2756 int *edgeptr;
2757 REAL *holeptr; /* unused */
2758 REAL xscale;
2759 REAL yscale;
2760 REAL xoffset;
2761 REAL yoffset;
2762 {
2763 int i;
2764 int index;
2765 REAL *point1, *point2;
2766
2767 index = dim;
2768 for (i = 1; i <= nodes; i++) {
2769 fprintf(polyfile, "%d %d %d 0 360 arc\nfill\n",
2770 (int) (nodeptr[index] * xscale + xoffset),
2771 (int) (nodeptr[index + 1] * yscale + yoffset),
2772 1 + (line_width >> 1));
2773 index += dim;
2774 }
2775 index = 2;
2776 for (i = 1; i <= edges; i++) {
2777 point1 = &nodeptr[edgeptr[index++] * dim];
2778 point2 = &nodeptr[edgeptr[index++] * dim];
2779 fprintf(polyfile, "%d %d moveto\n",
2780 (int) (point1[0] * xscale + xoffset),
2781 (int) (point1[1] * yscale + yoffset));
2782 fprintf(polyfile, "%d %d lineto\nstroke\n",
2783 (int) (point2[0] * xscale + xoffset),
2784 (int) (point2[1] * yscale + yoffset));
2785 }
2786 }
2787
print_ele(elefile,nodes,dim,elems,corners,nodeptr,eleptr,partition,shift,xscale,yscale,xoffset,yoffset,llcornerx,llcornery)2788 void print_ele(elefile, nodes, dim, elems, corners, nodeptr, eleptr,
2789 partition, shift,
2790 xscale, yscale, xoffset, yoffset, llcornerx, llcornery)
2791 FILE *elefile;
2792 int nodes; /* unused */
2793 int dim;
2794 int elems;
2795 int corners; /* unused */
2796 REAL *nodeptr;
2797 int *eleptr;
2798 int *partition;
2799 REAL *shift;
2800 REAL xscale;
2801 REAL yscale;
2802 REAL xoffset;
2803 REAL yoffset;
2804 int llcornerx;
2805 int llcornery;
2806 {
2807 int i, j;
2808 int index, colorindex;
2809 REAL shiftx, shifty;
2810 REAL *nowpoint;
2811
2812 index = 3;
2813 if ((partition != (int *) NULL) && !bw_ps) {
2814 fprintf(elefile, "0 0 0 setrgbcolor\n");
2815 fprintf(elefile, "%d %d moveto\n", llcornerx, llcornery);
2816 fprintf(elefile, "%d %d lineto\n", 612 - llcornerx, llcornery);
2817 fprintf(elefile, "%d %d lineto\n", 612 - llcornerx, 792 - llcornery);
2818 fprintf(elefile, "%d %d lineto\n", llcornerx, 792 - llcornery);
2819 fprintf(elefile, "fill\n");
2820 }
2821 for (i = 1; i <= elems; i++) {
2822 if ((partition != (int *) NULL) && !bw_ps) {
2823 colorindex = partition[i] & 63;
2824 fprintf(elefile, "%6.3f %6.3f %6.3f setrgbcolor\n",
2825 (REAL) rgb[colorindex].red / 65535.0,
2826 (REAL) rgb[colorindex].green / 65535.0,
2827 (REAL) rgb[colorindex].blue / 65535.0);
2828 }
2829 nowpoint = &nodeptr[eleptr[index + 2] * dim];
2830 if ((partition != (int *) NULL) && (explode || bw_ps)) {
2831 shiftx = shift[partition[i] << 1];
2832 shifty = shift[(partition[i] << 1) + 1];
2833 fprintf(elefile, "%d %d moveto\n",
2834 (int) ((nowpoint[0] + shiftx) * xscale + xoffset),
2835 (int) ((nowpoint[1] + shifty) * yscale + yoffset));
2836 for (j = 0; j < 3; j++) {
2837 nowpoint = &nodeptr[eleptr[index++] * dim];
2838 fprintf(elefile, "%d %d lineto\n",
2839 (int) ((nowpoint[0] + shiftx) * xscale + xoffset),
2840 (int) ((nowpoint[1] + shifty) * yscale + yoffset));
2841 }
2842 } else {
2843 fprintf(elefile, "%d %d moveto\n",
2844 (int) (nowpoint[0] * xscale + xoffset),
2845 (int) (nowpoint[1] * yscale + yoffset));
2846 for (j = 0; j < 3; j++) {
2847 nowpoint = &nodeptr[eleptr[index++] * dim];
2848 fprintf(elefile, "%d %d lineto\n",
2849 (int) (nowpoint[0] * xscale + xoffset),
2850 (int) (nowpoint[1] * yscale + yoffset));
2851 }
2852 }
2853 if (fillelem && !bw_ps) {
2854 fprintf(elefile, "gsave\nfill\ngrestore\n1 1 0 setrgbcolor\n");
2855 }
2856 fprintf(elefile, "stroke\n");
2857 }
2858 }
2859
print_edge(edgefile,nodes,dim,edges,nodeptr,edgeptr,normptr,xscale,yscale,xoffset,yoffset,llcornerx,llcornery)2860 void print_edge(edgefile, nodes, dim, edges, nodeptr, edgeptr, normptr,
2861 xscale, yscale, xoffset, yoffset, llcornerx, llcornery)
2862 FILE *edgefile;
2863 int nodes; /* unused */
2864 int dim;
2865 int edges;
2866 REAL *nodeptr;
2867 int *edgeptr;
2868 REAL *normptr;
2869 REAL xscale;
2870 REAL yscale;
2871 REAL xoffset;
2872 REAL yoffset;
2873 int llcornerx;
2874 int llcornery;
2875 {
2876 int i;
2877 int index;
2878 REAL *point1, *point2;
2879 REAL normx, normy;
2880 REAL normmult, normmultx, normmulty;
2881 REAL windowxmin, windowymin, windowxmax, windowymax;
2882
2883 index = 2;
2884 for (i = 1; i <= edges; i++) {
2885 point1 = &nodeptr[edgeptr[index++] * dim];
2886 if (edgeptr[index] == -1) {
2887 normx = normptr[index - 1];
2888 normy = normptr[index++];
2889 normmultx = 0.0;
2890 if (normx > 0) {
2891 windowxmax = ((REAL) (612 - llcornerx) - xoffset) / xscale;
2892 normmultx = (windowxmax - point1[0]) / normx;
2893 } else if (normx < 0) {
2894 windowxmin = ((REAL) llcornerx - xoffset) / xscale;
2895 normmultx = (windowxmin - point1[0]) / normx;
2896 }
2897 normmulty = 0.0;
2898 if (normy > 0) {
2899 windowymax = ((REAL) (792 - llcornery) - yoffset) / yscale;
2900 normmulty = (windowymax - point1[1]) / normy;
2901 } else if (normy < 0) {
2902 windowymin = ((REAL) llcornery - yoffset) / yscale;
2903 normmulty = (windowymin - point1[1]) / normy;
2904 }
2905 if (normmultx == 0.0) {
2906 normmult = normmulty;
2907 } else if (normmulty == 0.0) {
2908 normmult = normmultx;
2909 } else if (normmultx < normmulty) {
2910 normmult = normmultx;
2911 } else {
2912 normmult = normmulty;
2913 }
2914 if (normmult > 0.0) {
2915 fprintf(edgefile, "%d %d moveto\n",
2916 (int) (point1[0] * xscale + xoffset),
2917 (int) (point1[1] * yscale + yoffset));
2918 fprintf(edgefile, "%d %d lineto\nstroke\n",
2919 (int) ((point1[0] + normmult * normx) * xscale + xoffset),
2920 (int) ((point1[1] + normmult * normy) * yscale + yoffset));
2921 }
2922 } else {
2923 point2 = &nodeptr[edgeptr[index++] * dim];
2924 fprintf(edgefile, "%d %d moveto\n",
2925 (int) (point1[0] * xscale + xoffset),
2926 (int) (point1[1] * yscale + yoffset));
2927 fprintf(edgefile, "%d %d lineto\nstroke\n",
2928 (int) (point2[0] * xscale + xoffset),
2929 (int) (point2[1] * yscale + yoffset));
2930 }
2931 }
2932 }
2933
print_adj(adjfile,dim,subdomains,ptr,center,xscale,yscale,xoffset,yoffset,llcornerx,llcornery)2934 void print_adj(adjfile, dim, subdomains, ptr, center, xscale, yscale,
2935 xoffset, yoffset, llcornerx, llcornery)
2936 FILE *adjfile;
2937 int dim;
2938 int subdomains;
2939 int *ptr;
2940 REAL *center;
2941 REAL xscale;
2942 REAL yscale;
2943 REAL xoffset;
2944 REAL yoffset;
2945 int llcornerx;
2946 int llcornery;
2947 {
2948 int i, j;
2949 REAL *point1, *point2;
2950 int colorindex;
2951
2952 if (!bw_ps) {
2953 fprintf(adjfile, "0 0 0 setrgbcolor\n");
2954 fprintf(adjfile, "%d %d moveto\n", llcornerx, llcornery);
2955 fprintf(adjfile, "%d %d lineto\n", 612 - llcornerx, llcornery);
2956 fprintf(adjfile, "%d %d lineto\n", 612 - llcornerx, 792 - llcornery);
2957 fprintf(adjfile, "%d %d lineto\n", llcornerx, 792 - llcornery);
2958 fprintf(adjfile, "fill\n");
2959 fprintf(adjfile, "1 1 0 setrgbcolor\n");
2960 }
2961 for (i = 0; i < subdomains; i++) {
2962 for (j = i + 1; j < subdomains; j++) {
2963 if (ptr[i * subdomains + j]) {
2964 point1 = ¢er[i * dim];
2965 point2 = ¢er[j * dim];
2966 fprintf(adjfile, "%d %d moveto\n",
2967 (int) (point1[0] * xscale + xoffset),
2968 (int) (point1[1] * yscale + yoffset));
2969 fprintf(adjfile, "%d %d lineto\nstroke\n",
2970 (int) (point2[0] * xscale + xoffset),
2971 (int) (point2[1] * yscale + yoffset));
2972 }
2973 }
2974 }
2975 for (i = 0; i < subdomains; i++) {
2976 point1 = ¢er[i * dim];
2977 if (!bw_ps) {
2978 colorindex = i & 63;
2979 fprintf(adjfile, "%6.3f %6.3f %6.3f setrgbcolor\n",
2980 (REAL) rgb[colorindex].red / 65535.0,
2981 (REAL) rgb[colorindex].green / 65535.0,
2982 (REAL) rgb[colorindex].blue / 65535.0);
2983 fprintf(adjfile, "%d %d %d 0 360 arc\nfill\n",
2984 (int) (point1[0] * xscale + xoffset),
2985 (int) (point1[1] * yscale + yoffset),
2986 5 + (line_width >> 1));
2987 } else {
2988 fprintf(adjfile, "%d %d %d 0 360 arc\nfill\n",
2989 (int) (point1[0] * xscale + xoffset),
2990 (int) (point1[1] * yscale + yoffset),
2991 3 + (line_width >> 1));
2992 }
2993 }
2994 }
2995
print(inc,image,xmin,ymin,xmax,ymax,eps)2996 void print(inc, image, xmin, ymin, xmax, ymax, eps)
2997 int inc;
2998 int image;
2999 REAL xmin;
3000 REAL ymin;
3001 REAL xmax;
3002 REAL ymax;
3003 int eps;
3004 {
3005 REAL xxscale, yyscale, xxoffset, yyoffset;
3006 char psfilename[FILENAMESIZE];
3007 int llcornerx, llcornery;
3008 FILE *psfile;
3009
3010 if (image == NOTHING) {
3011 return;
3012 }
3013 if (!loaded[inc][image]) {
3014 return;
3015 }
3016 if ((image == PART) && (explode || bw_ps)) {
3017 xmin += (xmin - partcenter[inc][subdomains[inc] << 1]) * explosion;
3018 xmax += (xmax - partcenter[inc][subdomains[inc] << 1]) * explosion;
3019 ymin += (ymin - partcenter[inc][(subdomains[inc] << 1) + 1]) * explosion;
3020 ymax += (ymax - partcenter[inc][(subdomains[inc] << 1) + 1]) * explosion;
3021 }
3022 xxscale = (460.0 - (REAL) line_width) / (xmax - xmin);
3023 yyscale = (640.0 - (REAL) line_width) / (ymax - ymin);
3024 if (xxscale > yyscale) {
3025 xxscale = yyscale;
3026 llcornerx = (604 - (int) (yyscale * (xmax - xmin)) - line_width) >> 1;
3027 llcornery = 72;
3028 } else {
3029 yyscale = xxscale;
3030 llcornerx = 72;
3031 llcornery = (784 - (int) (xxscale * (ymax - ymin)) - line_width) >> 1;
3032 }
3033 xxoffset = 0.5 * (612.0 - xxscale * (xmax - xmin)) - xxscale * xmin +
3034 (line_width >> 1);
3035 yyoffset = 0.5 * (792.0 - yyscale * (ymax - ymin)) - yyscale * ymin +
3036 (line_width >> 1);
3037 switch(image) {
3038 case NODE:
3039 addps(nodefilename[inc], psfilename, eps);
3040 break;
3041 case POLY:
3042 addps(polyfilename[inc], psfilename, eps);
3043 break;
3044 case ELE:
3045 addps(elefilename[inc], psfilename, eps);
3046 break;
3047 case EDGE:
3048 addps(edgefilename[inc], psfilename, eps);
3049 break;
3050 case PART:
3051 addps(partfilename[inc], psfilename, eps);
3052 break;
3053 case ADJ:
3054 addps(adjfilename[inc], psfilename, eps);
3055 break;
3056 case VORO:
3057 addps(vedgefilename[inc], psfilename, eps);
3058 break;
3059 default:
3060 break;
3061 }
3062 if (print_head(psfilename, &psfile, llcornerx, llcornery, eps)) {
3063 return;
3064 }
3065 switch (image) {
3066 case NODE:
3067 print_node(psfile, nodes[inc], node_dim[inc], nodeptr[inc],
3068 xxscale, yyscale, xxoffset, yyoffset);
3069 break;
3070 case POLY:
3071 if (polynodes[inc] > 0) {
3072 print_poly(psfile, polynodes[inc], poly_dim[inc], polyedges[inc],
3073 polyholes[inc], polynodeptr[inc], polyedgeptr[inc],
3074 polyholeptr[inc], xxscale, yyscale, xxoffset, yyoffset);
3075 } else {
3076 print_poly(psfile, nodes[inc], node_dim[inc], polyedges[inc],
3077 polyholes[inc], nodeptr[inc], polyedgeptr[inc],
3078 polyholeptr[inc], xxscale, yyscale, xxoffset, yyoffset);
3079 }
3080 break;
3081 case ELE:
3082 print_ele(psfile, nodes[inc], node_dim[inc], elems[inc],
3083 ele_corners[inc], nodeptr[inc], eleptr[inc],
3084 (int *) NULL, (REAL *) NULL,
3085 xxscale, yyscale, xxoffset, yyoffset, llcornerx, llcornery);
3086 break;
3087 case EDGE:
3088 print_edge(psfile, nodes[inc], node_dim[inc], edges[inc],
3089 nodeptr[inc], edgeptr[inc], normptr[inc],
3090 xxscale, yyscale, xxoffset, yyoffset, llcornerx, llcornery);
3091 break;
3092 case PART:
3093 print_ele(psfile, nodes[inc], node_dim[inc], elems[inc],
3094 ele_corners[inc], nodeptr[inc], eleptr[inc],
3095 partpart[inc], partshift[inc],
3096 xxscale, yyscale, xxoffset, yyoffset, llcornerx, llcornery);
3097 break;
3098 case ADJ:
3099 print_adj(psfile, node_dim[inc], adjsubdomains[inc], adjptr[inc],
3100 partcenter[inc],
3101 xxscale, yyscale, xxoffset, yyoffset, llcornerx, llcornery);
3102 break;
3103 case VORO:
3104 print_edge(psfile, vnodes[inc], vnode_dim[inc], vedges[inc],
3105 vnodeptr[inc], vedgeptr[inc], vnormptr[inc],
3106 xxscale, yyscale, xxoffset, yyoffset, llcornerx, llcornery);
3107 break;
3108 default:
3109 break;
3110 }
3111 if (!eps) {
3112 fprintf(psfile, "showpage\n");
3113 }
3114 fclose(psfile);
3115 }
3116
main(argc,argv)3117 int main(argc, argv)
3118 int argc;
3119 char **argv;
3120 {
3121 REAL xmin = 0.0, ymin = 0.0, xmax = 0.0, ymax = 0.0;
3122 REAL xptr, yptr, xspan, yspan;
3123 int past_image;
3124 int new_image = 0;
3125 int new_inc = 0;
3126
3127 parsecommandline(argc, argv);
3128 showme_init();
3129 choose_image(start_inc, start_image);
3130 showme_window(argc, argv);
3131
3132 if (current_image != NOTHING) {
3133 xmin = xlo[current_inc][current_image];
3134 ymin = ylo[current_inc][current_image];
3135 xmax = xhi[current_inc][current_image];
3136 ymax = yhi[current_inc][current_image];
3137 zoom = 0;
3138 }
3139
3140 XMaskEvent(display, ExposureMask, &event);
3141 while (1) {
3142 switch (event.type) {
3143 case ButtonRelease:
3144 if (event.xany.window == quitwin) {
3145 XDestroyWindow(display, mainwindow);
3146 XCloseDisplay(display);
3147 return 0;
3148 } else if (event.xany.window == leftwin) {
3149 xspan = 0.25 * (xmax - xmin);
3150 xmin += xspan;
3151 xmax += xspan;
3152 draw(current_inc, current_image, xmin, ymin, xmax, ymax);
3153 } else if (event.xany.window == rightwin) {
3154 xspan = 0.25 * (xmax - xmin);
3155 xmin -= xspan;
3156 xmax -= xspan;
3157 draw(current_inc, current_image, xmin, ymin, xmax, ymax);
3158 } else if (event.xany.window == upwin) {
3159 yspan = 0.25 * (ymax - ymin);
3160 ymin -= yspan;
3161 ymax -= yspan;
3162 draw(current_inc, current_image, xmin, ymin, xmax, ymax);
3163 } else if (event.xany.window == downwin) {
3164 yspan = 0.25 * (ymax - ymin);
3165 ymin += yspan;
3166 ymax += yspan;
3167 draw(current_inc, current_image, xmin, ymin, xmax, ymax);
3168 } else if (event.xany.window == resetwin) {
3169 xmin = xlo[current_inc][current_image];
3170 ymin = ylo[current_inc][current_image];
3171 xmax = xhi[current_inc][current_image];
3172 ymax = yhi[current_inc][current_image];
3173 zoom = 0;
3174 draw(current_inc, current_image, xmin, ymin, xmax, ymax);
3175 } else if (event.xany.window == widthpluswin) {
3176 if (line_width < 100) {
3177 line_width++;
3178 XSetLineAttributes(display, linegc, line_width, LineSolid,
3179 CapRound, JoinRound);
3180 XSetLineAttributes(display, trianglegc, line_width, LineSolid,
3181 CapRound, JoinRound);
3182 draw(current_inc, current_image, xmin, ymin, xmax, ymax);
3183 }
3184 } else if (event.xany.window == widthminuswin) {
3185 if (line_width > 1) {
3186 line_width--;
3187 XSetLineAttributes(display, linegc, line_width, LineSolid,
3188 CapRound, JoinRound);
3189 XSetLineAttributes(display, trianglegc, line_width, LineSolid,
3190 CapRound, JoinRound);
3191 draw(current_inc, current_image, xmin, ymin, xmax, ymax);
3192 }
3193 } else if (event.xany.window == expwin) {
3194 if ((current_image == PART) && loaded[current_inc][PART]) {
3195 explode = !explode;
3196 draw(current_inc, current_image, xmin, ymin, xmax, ymax);
3197 }
3198 } else if (event.xany.window == exppluswin) {
3199 if ((current_image == PART) && loaded[current_inc][PART] && explode) {
3200 explosion += 0.125;
3201 findpartshift(subdomains[current_inc], explosion,
3202 partcenter[current_inc], partshift[current_inc]);
3203 draw(current_inc, current_image, xmin, ymin, xmax, ymax);
3204 }
3205 } else if (event.xany.window == expminuswin) {
3206 if ((current_image == PART) && loaded[current_inc][PART] && explode &&
3207 (explosion >= 0.125)) {
3208 explosion -= 0.125;
3209 findpartshift(subdomains[current_inc], explosion,
3210 partcenter[current_inc], partshift[current_inc]);
3211 draw(current_inc, current_image, xmin, ymin, xmax, ymax);
3212 }
3213 } else if (event.xany.window == fillwin) {
3214 if ((current_image == PART) && loaded[current_inc][PART]) {
3215 fillelem = !fillelem;
3216 draw(current_inc, current_image, xmin, ymin, xmax, ymax);
3217 }
3218 } else if (event.xany.window == pswin) {
3219 fill_button(pswin);
3220 XFlush(display);
3221 print(current_inc, current_image, xmin, ymin, xmax, ymax, 0);
3222 XClearWindow(display, pswin);
3223 XDrawString(display, pswin, fontgc, 2, 13, "PS", 2);
3224 } else if (event.xany.window == epswin) {
3225 fill_button(epswin);
3226 XFlush(display);
3227 print(current_inc, current_image, xmin, ymin, xmax, ymax, 1);
3228 XClearWindow(display, epswin);
3229 XDrawString(display, epswin, fontgc, 2, 13, "EPS", 3);
3230 } else if (event.xany.window == versionpluswin) {
3231 move_inc(1);
3232 loweriteration++;
3233 set_filenames(filename, loweriteration);
3234 if (current_inc == 1) {
3235 current_inc = 0;
3236 } else {
3237 current_image = NOTHING;
3238 XClearWindow(display, mainwindow);
3239 }
3240 draw_buttons();
3241 } else if (event.xany.window == versionminuswin) {
3242 if (loweriteration > 0) {
3243 move_inc(0);
3244 loweriteration--;
3245 set_filenames(filename, loweriteration);
3246 if (current_inc == 0) {
3247 current_inc = 1;
3248 } else {
3249 current_image = NOTHING;
3250 XClearWindow(display, mainwindow);
3251 }
3252 draw_buttons();
3253 }
3254 } else if ((event.xany.window == nodewin[0]) ||
3255 (event.xany.window == polywin[0]) ||
3256 (event.xany.window == elewin[0]) ||
3257 (event.xany.window == edgewin[0]) ||
3258 (event.xany.window == partwin[0]) ||
3259 (event.xany.window == adjwin[0]) ||
3260 (event.xany.window == voronoiwin[0]) ||
3261 (event.xany.window == nodewin[1]) ||
3262 (event.xany.window == polywin[1]) ||
3263 (event.xany.window == elewin[1]) ||
3264 (event.xany.window == edgewin[1]) ||
3265 (event.xany.window == partwin[1]) ||
3266 (event.xany.window == adjwin[1]) ||
3267 (event.xany.window == voronoiwin[1])) {
3268 if (event.xany.window == nodewin[0]) {
3269 new_inc = 0;
3270 new_image = NODE;
3271 }
3272 if (event.xany.window == polywin[0]) {
3273 new_inc = 0;
3274 new_image = POLY;
3275 }
3276 if (event.xany.window == elewin[0]) {
3277 new_inc = 0;
3278 new_image = ELE;
3279 }
3280 if (event.xany.window == edgewin[0]) {
3281 new_inc = 0;
3282 new_image = EDGE;
3283 }
3284 if (event.xany.window == partwin[0]) {
3285 new_inc = 0;
3286 new_image = PART;
3287 }
3288 if (event.xany.window == adjwin[0]) {
3289 new_inc = 0;
3290 new_image = ADJ;
3291 }
3292 if (event.xany.window == voronoiwin[0]) {
3293 new_inc = 0;
3294 new_image = VORO;
3295 }
3296 if (event.xany.window == nodewin[1]) {
3297 new_inc = 1;
3298 new_image = NODE;
3299 }
3300 if (event.xany.window == polywin[1]) {
3301 new_inc = 1;
3302 new_image = POLY;
3303 }
3304 if (event.xany.window == elewin[1]) {
3305 new_inc = 1;
3306 new_image = ELE;
3307 }
3308 if (event.xany.window == edgewin[1]) {
3309 new_inc = 1;
3310 new_image = EDGE;
3311 }
3312 if (event.xany.window == partwin[1]) {
3313 new_inc = 1;
3314 new_image = PART;
3315 }
3316 if (event.xany.window == adjwin[1]) {
3317 new_inc = 1;
3318 new_image = ADJ;
3319 }
3320 if (event.xany.window == voronoiwin[1]) {
3321 new_inc = 1;
3322 new_image = VORO;
3323 }
3324 past_image = current_image;
3325 if ((current_inc == new_inc) && (current_image == new_image)) {
3326 free_inc(new_inc);
3327 unload_inc(new_inc);
3328 }
3329 choose_image(new_inc, new_image);
3330 if ((past_image == NOTHING) && (current_image != NOTHING)) {
3331 xmin = xlo[current_inc][current_image];
3332 ymin = ylo[current_inc][current_image];
3333 xmax = xhi[current_inc][current_image];
3334 ymax = yhi[current_inc][current_image];
3335 zoom = 0;
3336 }
3337 draw(current_inc, current_image, xmin, ymin, xmax, ymax);
3338 } else {
3339 xptr = ((REAL) event.xbutton.x - xoffset) / xscale;
3340 yptr = ((REAL) event.xbutton.y - yoffset) / yscale;
3341 if ((current_image == PART) && loaded[current_inc][PART] && explode) {
3342 xptr = (xptr + partcenter[current_inc]
3343 [subdomains[current_inc] << 1]
3344 * explosion) / (1.0 + explosion);
3345 yptr = (yptr + partcenter[current_inc]
3346 [(subdomains[current_inc] << 1) + 1]
3347 * explosion) / (1.0 + explosion);
3348 }
3349 if ((event.xbutton.button == Button1)
3350 || (event.xbutton.button == Button3)) {
3351 if (event.xbutton.button == Button1) {
3352 xspan = 0.25 * (xmax - xmin);
3353 yspan = 0.25 * (ymax - ymin);
3354 zoom++;
3355 } else {
3356 xspan = xmax - xmin;
3357 yspan = ymax - ymin;
3358 zoom--;
3359 }
3360 xmin = xptr - xspan;
3361 ymin = yptr - yspan;
3362 xmax = xptr + xspan;
3363 ymax = yptr + yspan;
3364 draw(current_inc, current_image, xmin, ymin, xmax, ymax);
3365 } else if (event.xbutton.button == Button2) {
3366 printf("x = %.4g, y = %.4g\n", xptr, yptr);
3367 }
3368 }
3369 break;
3370 case DestroyNotify:
3371 XDestroyWindow(display, mainwindow);
3372 XCloseDisplay(display);
3373 return 0;
3374 case ConfigureNotify:
3375 if ((width != event.xconfigure.width) ||
3376 (height != event.xconfigure.height - PANELHEIGHT)) {
3377 width = event.xconfigure.width;
3378 height = event.xconfigure.height - PANELHEIGHT;
3379 draw(current_inc, current_image, xmin, ymin, xmax, ymax);
3380 while (XCheckMaskEvent(display, ExposureMask, &event));
3381 }
3382 break;
3383 case Expose:
3384 draw(current_inc, current_image, xmin, ymin, xmax, ymax);
3385 while (XCheckMaskEvent(display, ExposureMask, &event));
3386 break;
3387 default:
3388 break;
3389 }
3390 XNextEvent(display, &event);
3391 }
3392 }
3393