ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/gen/genprism.c
(Generate patch)

Comparing ray/src/gen/genprism.c (file contents):
Revision 2.4 by greg, Mon Aug 2 14:23:04 1993 UTC vs.
Revision 2.13 by greg, Sat Jul 25 19:18:01 2020 UTC

# Line 1 | Line 1
1 /* Copyright (c) 1987 Regents of the University of California */
2
1   #ifndef lint
2 < static char SCCSid[] = "$SunId$ LBL";
2 > static const char       RCSid[] = "$Id$";
3   #endif
6
4   /*
5   *  genprism.c - generate a prism.
6   *              2D vertices in the xy plane are given on the
7   *              command line or from a file.  Their order together
8   *              with the extrude direction will determine surface
9   *              orientation.
13 *
14 *     8/24/87
10   */
11  
12 < #include  <stdio.h>
13 <
12 > #include  "rtio.h"
13 > #include  <stdlib.h>
14   #include  <math.h>
20
15   #include  <ctype.h>
16  
17   #define  MAXVERT        1024            /* maximum # vertices */
18  
19 < #ifdef  DCL_ATOF
26 < extern double  atof();
27 < #endif
19 > #define  FTINY          1e-6
20  
21   char  *pmtype;          /* material type */
22   char  *pname;           /* name */
23  
24   double  lvect[3] = {0.0, 0.0, 1.0};
25 + int     lvdir = 1;
26 + double  llen = 1.0;
27  
28   double  vert[MAXVERT][2];
29   int  nverts = 0;
30  
31 + double  u[MAXVERT][2];          /* edge unit vectors */
32 + double  a[MAXVERT];             /* corner trim sizes */
33 +
34   int  do_ends = 1;               /* include end caps */
35   int  iscomplete = 0;            /* polygon is already completed */
36 + double  crad = 0.0;             /* radius for corners (sign from lvdir) */
37  
38 + static double  compute_rounding(void);
39  
40 < main(argc, argv)
41 < int  argc;
42 < char  **argv;
40 >
41 > static void
42 > readverts(char *fname)          /* read vertices from a file */
43   {
44 +        FILE  *fp;
45 +
46 +        if (fname == NULL)
47 +                fp = stdin;
48 +        else if ((fp = fopen(fname, "r")) == NULL) {
49 +                fprintf(stderr, "%s: cannot open\n", fname);
50 +                exit(1);
51 +        }
52 +        while (fscanf(fp, "%lf %lf", &vert[nverts][0], &vert[nverts][1]) == 2)
53 +                nverts++;
54 +        fclose(fp);
55 + }
56 +
57 +
58 + static void
59 + side(int n0, int n1)                    /* print single side */
60 + {
61 +        printf("\n%s polygon %s.%d\n", pmtype, pname, n0+1);
62 +        printf("0\n0\n12\n");
63 +        printf("\t%18.12g\t%18.12g\t%18.12g\n", vert[n0][0],
64 +                        vert[n0][1], 0.0);
65 +        printf("\t%18.12g\t%18.12g\t%18.12g\n", vert[n0][0]+lvect[0],
66 +                        vert[n0][1]+lvect[1], lvect[2]);
67 +        printf("\t%18.12g\t%18.12g\t%18.12g\n", vert[n1][0]+lvect[0],
68 +                        vert[n1][1]+lvect[1], lvect[2]);
69 +        printf("\t%18.12g\t%18.12g\t%18.12g\n", vert[n1][0],
70 +                        vert[n1][1], 0.0);
71 + }
72 +
73 +
74 + static void
75 + rside(int n0, int n1)                   /* print side with rounded edge */
76 + {
77 +        double  s, c, t[3];
78 +
79 +                                        /* compute tanget offset vector */
80 +        s = lvdir*(lvect[1]*u[n1][0] - lvect[0]*u[n1][1])/llen;
81 +        if (s < -FTINY || s > FTINY) {
82 +                c = sqrt(1. - s*s);
83 +                t[0] = (c - 1.)*u[n1][1];
84 +                t[1] = (1. - c)*u[n1][0];
85 +                t[2] = s;
86 +        } else
87 +                t[0] = t[1] = t[2] = 0.;
88 +                                        /* output side */
89 +        printf("\n%s polygon %s.%d\n", pmtype, pname, n0+1);
90 +        printf("0\n0\n12\n");
91 +        printf("\t%18.12g\t%18.12g\t%18.12g\n",
92 +                        vert[n0][0] + crad*(t[0] - a[n0]*u[n1][0]),
93 +                        vert[n0][1] + crad*(t[1] - a[n0]*u[n1][1]),
94 +                        crad*(t[2] + 1.));
95 +        printf("\t%18.12g\t%18.12g\t%18.12g\n",
96 +                        vert[n0][0] + lvect[0] + crad*(t[0] - a[n0]*u[n1][0]),
97 +                        vert[n0][1] + lvect[1] + crad*(t[1] - a[n0]*u[n1][1]),
98 +                        lvect[2] + crad*(t[2] - 1.));
99 +        printf("\t%18.12g\t%18.12g\t%18.12g\n",
100 +                        vert[n1][0] + lvect[0] + crad*(t[0] + a[n1]*u[n1][0]),
101 +                        vert[n1][1] + lvect[1] + crad*(t[1] + a[n1]*u[n1][1]),
102 +                        lvect[2] + crad*(t[2] - 1.));
103 +        printf("\t%18.12g\t%18.12g\t%18.12g\n",
104 +                        vert[n1][0] + crad*(t[0] + a[n1]*u[n1][0]),
105 +                        vert[n1][1] + crad*(t[1] + a[n1]*u[n1][1]),
106 +                        crad*(t[2] + 1.));
107 +                                        /* output joining edge */
108 +        if (lvdir*a[n1] < 0.)
109 +                return;
110 +        printf("\n%s cylinder %s.e%d\n", pmtype, pname, n0+1);
111 +        printf("0\n0\n7\n");
112 +        printf("\t%18.12g\t%18.12g\t%18.12g\n",
113 +                vert[n1][0] + crad*(a[n1]*u[n1][0] - u[n1][1]),
114 +                vert[n1][1] + crad*(a[n1]*u[n1][1] + u[n1][0]),
115 +                crad);
116 +        printf("\t%18.12g\t%18.12g\t%18.12g\n",
117 +                vert[n1][0] + lvect[0] + crad*(a[n1]*u[n1][0] - u[n1][1]),
118 +                vert[n1][1] + lvect[1] + crad*(a[n1]*u[n1][1] + u[n1][0]),
119 +                lvect[2] - crad);
120 +        printf("\t%18.12g\n", lvdir*crad);
121 + }
122 +
123 +
124 + static double
125 + compute_rounding(void)          /* compute vectors for rounding operations */
126 + {
127 +        int  i;
128 +        double  *v0, *v1;
129 +        double  l, asum;
130 +
131 +        v0 = vert[nverts-1];
132 +        for (i = 0; i < nverts; i++) {          /* compute u[*] */
133 +                v1 = vert[i];
134 +                u[i][0] = v0[0] - v1[0];
135 +                u[i][1] = v0[1] - v1[1];
136 +                l = sqrt(u[i][0]*u[i][0] + u[i][1]*u[i][1]);
137 +                if (l <= FTINY) {
138 +                        fprintf(stderr, "Degenerate side in prism\n");
139 +                        exit(1);
140 +                }
141 +                u[i][0] /= l;
142 +                u[i][1] /= l;
143 +                v0 = v1;
144 +        }
145 +        asum = 0.;
146 +        v1 = u[0];
147 +        for (i = nverts; i--; ) {               /* compute a[*] */
148 +                v0 = u[i];
149 +                l = v0[0]*v1[0] + v0[1]*v1[1];
150 +                if (1+l <= FTINY) {
151 +                        fprintf(stderr, "Overlapping sides in prism\n");
152 +                        exit(1);
153 +                }
154 +                if (1-l <= 0.)
155 +                        a[i] = 0.;
156 +                else {
157 +                        a[i] = sqrt((1-l)/(1+l));
158 +                        asum += l = v1[0]*v0[1]-v1[1]*v0[0];
159 +                        if (l < 0.)
160 +                                a[i] = -a[i];
161 +                }
162 +                v1 = v0;
163 +        }
164 +        return(asum*.5);
165 + }
166 +
167 +
168 + static void
169 + printends(void)                 /* print ends of prism */
170 + {
171 +        int  i;
172 +                                                /* bottom face */
173 +        printf("\n%s polygon %s.b\n", pmtype, pname);
174 +        printf("0\n0\n%d\n", nverts*3);
175 +        for (i = 0; i < nverts; i++)
176 +                printf("\t%18.12g\t%18.12g\t%18.12g\n", vert[i][0],
177 +                                vert[i][1], 0.0);
178 +                                                /* top face */
179 +        printf("\n%s polygon %s.t\n", pmtype, pname);
180 +        printf("0\n0\n%d\n", nverts*3);
181 +        for (i = nverts; i--; )
182 +                printf("\t%18.12g\t%18.12g\t%18.12g\n", vert[i][0]+lvect[0],
183 +                                vert[i][1]+lvect[1], lvect[2]);
184 + }
185 +
186 +
187 + static void
188 + printrends(void)                /* print ends of prism with rounding */
189 + {
190 +        int  i;
191 +        double  c0[3], c1[3], cl[3];
192 +                                                /* bottom face */
193 +        printf("\n%s polygon %s.b\n", pmtype, pname);
194 +        printf("0\n0\n%d\n", nverts*3);
195 +        for (i = 0; i < nverts; i++) {
196 +                printf("\t%18.12g\t%18.12g\t%18.12g\n",
197 +                        vert[i][0] + crad*(a[i]*u[i][0] - u[i][1]),
198 +                        vert[i][1] + crad*(a[i]*u[i][1] + u[i][0]),
199 +                        0.0);
200 +        }
201 +                                                /* top face */
202 +        printf("\n%s polygon %s.t\n", pmtype, pname);
203 +        printf("0\n0\n%d\n", nverts*3);
204 +        for (i = nverts; i--; ) {
205 +                printf("\t%18.12g\t%18.12g\t%18.12g\n",
206 +                vert[i][0] + lvect[0] + crad*(a[i]*u[i][0] - u[i][1]),
207 +                vert[i][1] + lvect[1] + crad*(a[i]*u[i][1] + u[i][0]),
208 +                lvect[2]);
209 +        }
210 +                                                /* bottom corners and edges */
211 +        c0[0] = cl[0] = vert[nverts-1][0] +
212 +                        crad*(a[nverts-1]*u[nverts-1][0] - u[nverts-1][1]);
213 +        c0[1] = cl[1] = vert[nverts-1][1] +
214 +                        crad*(a[nverts-1]*u[nverts-1][1] + u[nverts-1][0]);
215 +        c0[2] = cl[2] = crad;
216 +        for (i = 0; i < nverts; i++) {
217 +                if (i < nverts-1) {
218 +                        c1[0] = vert[i][0] + crad*(a[i]*u[i][0] - u[i][1]);
219 +                        c1[1] = vert[i][1] + crad*(a[i]*u[i][1] + u[i][0]);
220 +                        c1[2] = crad;
221 +                } else {
222 +                        c1[0] = cl[0]; c1[1] = cl[1]; c1[2] = cl[2];
223 +                }
224 +                if (lvdir*a[i] > 0.) {
225 +                        printf("\n%s sphere %s.bc%d\n", pmtype, pname, i+1);
226 +                        printf("0\n0\n4 %18.12g %18.12g %18.12g %18.12g\n",
227 +                                c1[0], c1[1], c1[2], lvdir*crad);
228 +                }
229 +                printf("\n%s cylinder %s.be%d\n", pmtype, pname, i+1);
230 +                printf("0\n0\n7\n");
231 +                printf("\t%18.12g\t%18.12g\t%18.12g\n", c0[0], c0[1], c0[2]);
232 +                printf("\t%18.12g\t%18.12g\t%18.12g\n", c1[0], c1[1], c1[2]);
233 +                printf("\t%18.12g\n", lvdir*crad);
234 +                c0[0] = c1[0]; c0[1] = c1[1]; c0[2] = c1[2];
235 +        }
236 +                                                /* top corners and edges */
237 +        c0[0] = cl[0] = vert[nverts-1][0] + lvect[0] +
238 +                        crad*(a[nverts-1]*u[nverts-1][0] - u[nverts-1][1]);
239 +        c0[1] = cl[1] = vert[nverts-1][1] + lvect[1] +
240 +                        crad*(a[nverts-1]*u[nverts-1][1] + u[nverts-1][0]);
241 +        c0[2] = cl[2] = lvect[2] - crad;
242 +        for (i = 0; i < nverts; i++) {
243 +                if (i < nverts-1) {
244 +                        c1[0] = vert[i][0] + lvect[0] +
245 +                                        crad*(a[i]*u[i][0] - u[i][1]);
246 +                        c1[1] = vert[i][1] + lvect[1] +
247 +                                        crad*(a[i]*u[i][1] + u[i][0]);
248 +                        c1[2] = lvect[2] - crad;
249 +                } else {
250 +                        c1[0] = cl[0]; c1[1] = cl[1]; c1[2] = cl[2];
251 +                }
252 +                if (lvdir*a[i] > 0.) {
253 +                        printf("\n%s sphere %s.tc%d\n", pmtype, pname, i+1);
254 +                        printf("0\n0\n4 %18.12g %18.12g %18.12g %18.12g\n",
255 +                                c1[0], c1[1], c1[2], lvdir*crad);
256 +                }
257 +                printf("\n%s cylinder %s.te%d\n", pmtype, pname, i+1);
258 +                printf("0\n0\n7\n");
259 +                printf("\t%18.12g\t%18.12g\t%18.12g\n", c0[0], c0[1], c0[2]);
260 +                printf("\t%18.12g\t%18.12g\t%18.12g\n", c1[0], c1[1], c1[2]);
261 +                printf("\t%18.12g\n", lvdir*crad);
262 +                c0[0] = c1[0]; c0[1] = c1[1]; c0[2] = c1[2];
263 +        }
264 + }
265 +
266 +
267 + static void
268 + printsides(int round)           /* print prism sides */
269 + {
270 +        int  i;
271 +
272 +        for (i = 0; i < nverts-1; i++)
273 +                if (round)
274 +                        rside(i, i+1);
275 +                else
276 +                        side(i, i+1);
277 +        if (!iscomplete) {
278 +                if (round)
279 +                        rside(nverts-1, 0);
280 +                else
281 +                        side(nverts-1, 0);
282 +        }
283 + }
284 +
285 +
286 + int
287 + main(int argc, char **argv)
288 + {
289          int  an;
290          
291          if (argc < 4)
# Line 79 | Line 323 | char  **argv;
323                          lvect[0] = atof(argv[++an]);
324                          lvect[1] = atof(argv[++an]);
325                          lvect[2] = atof(argv[++an]);
326 +                        if (lvect[2] < -FTINY)
327 +                                lvdir = -1;
328 +                        else if (lvect[2] > FTINY)
329 +                                lvdir = 1;
330 +                        else {
331 +                                fprintf(stderr,
332 +                                        "%s: illegal extrusion vector\n",
333 +                                                argv[0]);
334 +                                exit(1);
335 +                        }
336 +                        llen = sqrt(lvect[0]*lvect[0] + lvect[1]*lvect[1] +
337 +                                        lvect[2]*lvect[2]);
338                          break;
339 +                case 'r':                               /* radius */
340 +                        crad = atof(argv[++an]);
341 +                        break;
342                  case 'e':                               /* ends */
343                          do_ends = !do_ends;
344                          break;
# Line 90 | Line 349 | char  **argv;
349                          goto userr;
350                  }
351          }
352 <
353 <        printhead(argc, argv);
354 <
355 <        if (do_ends)
356 <                printends();
357 <        printsides();
358 <
352 >        if (crad > FTINY) {
353 >                if (crad > lvdir*lvect[2]) {
354 >                        fprintf(stderr, "%s: rounding greater than height\n",
355 >                                        argv[0]);
356 >                        exit(1);
357 >                }
358 >                crad *= lvdir;          /* simplifies formulas */
359 >                compute_rounding();
360 >                fputs("# ", stdout);
361 >                printargs(argc, argv, stdout);
362 >                if (do_ends)
363 >                        printrends();
364 >                printsides(1);
365 >        } else {
366 >                fputs("# ", stdout);
367 >                printargs(argc, argv, stdout);
368 >                if (do_ends)
369 >                        printends();
370 >                printsides(0);
371 >        }
372          return(0);
373   userr:
374          fprintf(stderr, "Usage: %s material name ", argv[0]);
375          fprintf(stderr, "{ - | vfile | N v1 v2 .. vN } ");
376 <        fprintf(stderr, "[-l lvect][-c][-e]\n");
377 <        exit(1);
376 >        fprintf(stderr, "[-l lvect][-r radius][-c][-e]\n");
377 >        return(1);
378   }
379  
108
109 readverts(fname)                /* read vertices from a file */
110 char  *fname;
111 {
112        FILE  *fp;
113
114        if (fname == NULL)
115                fp = stdin;
116        else if ((fp = fopen(fname, "r")) == NULL) {
117                fprintf(stderr, "%s: cannot open\n", fname);
118                exit(1);
119        }
120        while (fscanf(fp, "%lf %lf", &vert[nverts][0], &vert[nverts][1]) == 2)
121                nverts++;
122        fclose(fp);
123 }
124
125
126 printends()                     /* print ends of prism */
127 {
128        register int  i;
129
130        printf("\n%s polygon %s.b\n", pmtype, pname);
131        printf("0\n0\n%d\n", nverts*3);
132        for (i = 0; i < nverts; i++) {
133                printf("\t%18.12g\t%18.12g\t%18.12g\n",
134                                vert[i][0],
135                                vert[i][1],
136                                0.0);
137        }
138        printf("\n%s polygon %s.e\n", pmtype, pname);
139        printf("0\n0\n%d\n", nverts*3);
140        for (i = nverts-1; i >= 0; i--) {
141                printf("\t%18.12g\t%18.12g\t%18.12g\n",
142                                vert[i][0]+lvect[0],
143                                vert[i][1]+lvect[1],
144                                lvect[2]);
145        }
146 }
147
148
149 printsides()                    /* print prism sides */
150 {
151        register int  i;
152
153        for (i = 0; i < nverts-1; i++)
154                side(i, i+1);
155        if (!iscomplete)
156                side(nverts-1, 0);
157 }
158
159
160 side(n1, n2)                    /* print single side */
161 register int  n1, n2;
162 {
163        printf("\n%s polygon %s.%d\n", pmtype, pname, n1+1);
164        printf("0\n0\n12\n");
165        printf("\t%18.12g\t%18.12g\t%18.12g\n",
166                        vert[n1][0],
167                        vert[n1][1],
168                        0.0);
169        printf("\t%18.12g\t%18.12g\t%18.12g\n",
170                        vert[n1][0]+lvect[0],
171                        vert[n1][1]+lvect[1],
172                        lvect[2]);
173        printf("\t%18.12g\t%18.12g\t%18.12g\n",
174                        vert[n2][0]+lvect[0],
175                        vert[n2][1]+lvect[1],
176                        lvect[2]);
177        printf("\t%18.12g\t%18.12g\t%18.12g\n",
178                        vert[n2][0],
179                        vert[n2][1],
180                        0.0);
181 }
182
183
184 printhead(ac, av)               /* print command header */
185 register int  ac;
186 register char  **av;
187 {
188        putchar('#');
189        while (ac--) {
190                putchar(' ');
191                fputs(*av++, stdout);
192        }
193        putchar('\n');
194 }

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines