ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/util/contour.c
Revision: 1.4
Committed: Fri Mar 26 23:34:23 2004 UTC (20 years ago) by schorsch
Content type: text/plain
Branch: MAIN
CVS Tags: rad3R6P1, rad3R6
Changes since 1.3: +52 -28 lines
Log Message:
Continued ANSIfication.

File Contents

# Content
1 #ifndef lint
2 static const char RCSid[] = "$Id: contour.c,v 1.3 2003/02/22 02:07:30 greg Exp $";
3 #endif
4 /*
5 * contour.c - program to make contour plots, mappings from 3 to
6 * 2 dimenstions.
7 *
8 * 8/22/86
9 */
10
11 #include <stdio.h>
12
13 #include <stdlib.h>
14
15 #include <ctype.h>
16
17 #define MAXPTS 2048 /* maximum number of input points */
18
19 typedef float DATAPT[3];
20
21 DATAPT xyz[MAXPTS]; /* the input data */
22 int xyzsiz;
23
24 double zmin = 1e20; /* z minimum */
25 double zmax = -1e20; /* z maximum */
26 int minset = 0; /* user set minimum? */
27 int maxset = 0; /* user set maximum? */
28
29 int ncurves = 6; /* number of contours */
30
31 static void getdata(FILE *fp);
32 static int xycmp(const void* p1, const void* p2);
33 static void contour(int n);
34 static void crossings(double z);
35 static void cross(DATAPT p0, DATAPT p1, double z);
36
37
38 int
39 main(
40 int argc,
41 char *argv[]
42 )
43 {
44 FILE *fp;
45 int i;
46
47 for (i = 1; i < argc && argv[i][0] == '-'; i++)
48 switch (argv[i][1]) {
49 case 'n':
50 ncurves = atoi(argv[++i]);
51 break;
52 case 'l':
53 zmin = atof(argv[++i]);
54 minset = 1;
55 break;
56 case 'u':
57 zmax = atof(argv[++i]);
58 maxset = 1;
59 break;
60 default:
61 fprintf(stderr, "%s: unknown option: %s\n",
62 argv[0], argv[i]);
63 exit(1);
64 }
65
66 if (i == argc)
67 getdata(stdin);
68 else if (i < argc-1) {
69 fprintf(stderr, "%s: too many file names\n", argv[0]);
70 exit(1);
71 } else if ((fp = fopen(argv[i], "r")) == NULL) {
72 fprintf(stderr, "%s: file not found\n", argv[i]);
73 exit(1);
74 } else
75 getdata(fp);
76
77 qsort(xyz, xyzsiz, sizeof(DATAPT), xycmp); /* sort data */
78
79 for (i = 0; i < ncurves; i++) /* do contours */
80 contour(i);
81
82 exit(0);
83 }
84
85
86 static void
87 getdata( /* read input data */
88 FILE *fp
89 )
90 {
91 register int i;
92
93 while (xyzsiz < MAXPTS) {
94 for (i = 0; i < 3; i++)
95 if (fscanf(fp, "%f", &xyz[xyzsiz][i]) != 1)
96 break;
97 if (i != 3)
98 break;
99 if (xyz[xyzsiz][2] < zmin) {
100 if (minset) {
101 continue;
102 } else {
103 zmin = xyz[xyzsiz][2];
104 }
105 }
106 if (xyz[xyzsiz][2] > zmax) {
107 if (maxset) {
108 continue;
109 } else {
110 zmax = xyz[xyzsiz][2];
111 }
112 }
113 xyzsiz++;
114 }
115 if (!feof(fp)) {
116 fprintf(stderr, "Error reading input data\n");
117 exit(1);
118 }
119 }
120
121
122 static int
123 xycmp( /* -1 if p1 < p2, 0 if equal, 1 if p1 > p2 */
124 const void* p1,
125 const void* p2
126 )
127 {
128 const float* pp1 = p1;
129 const float* pp2 = p2;
130 if (pp1[0] > pp2[0])
131 return(1);
132 if (pp1[0] < pp2[0])
133 return(-1);
134 if (pp1[1] > pp2[1])
135 return(1);
136 if (pp1[1] < pp2[1])
137 return(-1);
138 return(0);
139 }
140
141
142 static void
143 contour( /* make contour n */
144 int n
145 )
146 {
147 double z;
148
149 z = (n+.5)*(zmax-zmin)/ncurves + zmin;
150 printf("%clabel=\"%f\";\n", n+'A', z);
151 printf("%cdata=\n", n+'A');
152 crossings(z);
153 printf(";\n");
154 }
155
156
157 static void
158 crossings( /* find crossings for z */
159 double z
160 )
161 {
162 register DATAPT *p0, *p1;
163
164 p0 = p1 = xyz;
165 while (p0 < xyz+xyzsiz-2) {
166 while (p1 < xyz+xyzsiz-2) {
167 if (p0[0][0] < p1[0][0] && p0[0][1] <= p1[0][1])
168 break;
169 p1++;
170 }
171 if (p0[0][0] == p0[1][0])
172 cross(p0[0], p0[1], z);
173 if (p0[0][1] == p1[0][1])
174 cross(p0[0], p1[0], z);
175 p0++;
176 }
177 }
178
179
180 static void
181 cross( /* mark crossing between p0 and p1 */
182 register DATAPT p0,
183 register DATAPT p1,
184 double z
185 )
186 {
187 if (p1[2] - p0[2] == 0.0)
188 return;
189 z = (z - p0[2])/(p1[2] - p0[2]);
190 if (z < 0.0 || z >= 1.0)
191 return;
192 printf("%e\t", p0[0]*(1.0 - z) + p1[0]*z);
193 printf("%e\n", p0[1]*(1.0 - z) + p1[1]*z);
194 }