ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/util/eplus_adduvf.c
Revision: 2.14
Committed: Thu Feb 27 22:42:53 2014 UTC (10 years, 1 month ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 2.13: +7 -7 lines
Log Message:
Bug fix in last change for subsurfaces

File Contents

# User Rev Content
1 greg 2.1 #ifndef lint
2 greg 2.14 static const char RCSid[] = "$Id: eplus_adduvf.c,v 2.13 2014/02/27 20:39:04 greg Exp $";
3 greg 2.1 #endif
4     /*
5     * Add User View Factors to EnergyPlus Input Data File
6     *
7     * G.Ward for LBNL
8     */
9    
10     #include <stdlib.h>
11 greg 2.3 #include "rtio.h"
12     #include "rtmath.h"
13     #include "random.h"
14 greg 2.1 #include "eplus_idf.h"
15     #include "triangulate.h"
16 greg 2.2 #include "rtprocess.h"
17 greg 2.1
18     #ifndef NSAMPLES
19 greg 2.7 #define NSAMPLES 80000 /* default number of samples */
20 greg 2.1 #endif
21    
22 greg 2.11 #define SURF_EPS 0.0005 /* surface testing epsilon */
23    
24 greg 2.2 char *progname; /* global argv[0] */
25    
26 greg 2.7 int nsamps = NSAMPLES; /* number of samples to use */
27    
28 greg 2.5 char temp_octree[128]; /* temporary octree */
29 greg 2.2
30     const char UVF_PNAME[] =
31 greg 2.10 "ZoneProperty:UserViewFactors:bySurfaceName";
32 greg 2.1
33 greg 2.11 const char SUBSURF_PNAME[] =
34     "FenestrationSurface:Detailed";
35    
36 greg 2.1 const char ADD_HEADER[] =
37 greg 2.5 "\n!+++ User View Factors computed by Radiance +++!\n\n";
38 greg 2.1
39     #define NAME_FLD 1 /* name field always first? */
40    
41 greg 2.11 #define SS_BASE_FLD 4 /* subsurface base surface */
42     #define SS_VERT_FLD 10 /* subsurface vertex count */
43    
44 greg 2.1 typedef struct {
45 greg 2.9 const char *pname; /* object type name */
46 greg 2.1 short zone_fld; /* zone field index */
47     short vert_fld; /* vertex field index */
48     } SURF_PTYPE; /* surface type we're interested in */
49    
50     const SURF_PTYPE surf_type[] = {
51     {"BuildingSurface:Detailed", 4, 10},
52 greg 2.5 {"Floor:Detailed", 3, 9},
53     {"RoofCeiling:Detailed", 3, 9},
54     {"Wall:Detailed", 3, 9},
55 greg 2.1 {NULL}
56 greg 2.7 }; /* IDF surface types */
57 greg 2.1
58 greg 2.2 typedef struct s_zone {
59 greg 2.1 const char *zname; /* zone name */
60     struct s_zone *next; /* next zone in list */
61 greg 2.2 int nsurf; /* surface count */
62 greg 2.11 int ntotal; /* surfaces+subsurfaces */
63 greg 2.9 IDF_OBJECT *pfirst; /* first matching object */
64 greg 2.11 IDF_OBJECT *plast; /* last before subsurfaces */
65 greg 2.12 float *area_redu; /* subsurface area per surf. */
66 greg 2.2 } ZONE; /* a list of collected zone surfaces */
67    
68     ZONE *zone_list = NULL; /* our list of zones */
69 greg 2.1
70 greg 2.11 LUTAB zonesurf_lut = /* zone surface lookup table */
71     LU_SINIT(NULL,NULL);
72    
73 greg 2.1 IDF_LOADED *our_idf = NULL; /* loaded/modified IDF */
74    
75 greg 2.7 typedef struct {
76 greg 2.11 FVECT norm; /* surface normal */
77     double area; /* surface area */
78     int nv; /* number of vertices */
79     FVECT vl[3]; /* vertex list (extends struct) */
80     } SURFACE; /* a polygonal surface */
81    
82     typedef struct {
83 greg 2.7 FVECT sdir[3]; /* UVW unit sampling vectors */
84     double poff; /* W-offset for plane of polygon */
85     double area_left; /* area left to sample */
86     int samp_left; /* remaining samples */
87     int wd; /* output file descriptor */
88     } POLYSAMP; /* structure for polygon sampling */
89    
90 greg 2.1 /* Create a new zone and push to top of our list */
91 greg 2.2 static ZONE *
92 greg 2.9 new_zone(const char *zname, IDF_OBJECT *param)
93 greg 2.1 {
94 greg 2.2 ZONE *znew = (ZONE *)malloc(sizeof(ZONE));
95 greg 2.1
96     if (znew == NULL)
97     return(NULL);
98     znew->zname = zname; /* assumes static */
99     znew->next = zone_list;
100     znew->pfirst = znew->plast = param;
101 greg 2.11 znew->ntotal = znew->nsurf = 1;
102 greg 2.12 znew->area_redu = NULL;
103 greg 2.1 return(zone_list = znew);
104     }
105    
106     /* Add the detailed surface (polygon) to the named zone */
107 greg 2.2 static ZONE *
108 greg 2.9 add2zone(IDF_OBJECT *param, const char *zname)
109 greg 2.1 {
110 greg 2.11 IDF_FIELD *nfp = idf_getfield(param, NAME_FLD);
111     ZONE *zptr;
112     LUENT *lep;
113 greg 2.1
114 greg 2.11 if (nfp == NULL) {
115     fputs(progname, stderr);
116     fputs(": surface missing name field!\n", stderr);
117     return(NULL);
118     }
119 greg 2.1 for (zptr = zone_list; zptr != NULL; zptr = zptr->next)
120     if (!strcmp(zptr->zname, zname))
121     break;
122 greg 2.11 if (zptr == NULL) {
123     zptr = new_zone(zname, param);
124     } else { /* keep surfaces together */
125     if (!idf_movobject(our_idf, param, zptr->plast))
126     return(NULL);
127     zptr->plast = param;
128     zptr->nsurf++;
129     zptr->ntotal++;
130     }
131     /* add to lookup table */
132     lep = lu_find(&zonesurf_lut, nfp->val);
133     if (lep == NULL)
134     return(NULL);
135     if (lep->data != NULL) {
136     fputs(progname, stderr);
137     fputs(": duplicate surface name '", stderr);
138     fputs(nfp->val, stderr);
139     fputs("'\n", stderr);
140     return(NULL);
141     }
142     lep->key = nfp->val;
143     lep->data = (char *)zptr;
144     return(zptr);
145     }
146    
147     /* Add a subsurface by finding its base surface and the corresponding zone */
148     static ZONE *
149     add_subsurf(IDF_OBJECT *param)
150     {
151     IDF_FIELD *bfp = idf_getfield(param, SS_BASE_FLD);
152     ZONE *zptr;
153     LUENT *lep;
154    
155     if (bfp == NULL) {
156     fputs(progname, stderr);
157     fputs(": missing base field name in subsurface!\n", stderr);
158     return(NULL);
159     }
160     lep = lu_find(&zonesurf_lut, bfp->val); /* find base zone */
161     if (lep == NULL || lep->data == NULL) {
162     fputs(progname, stderr);
163     fputs(": cannot find referenced base surface '", stderr);
164     fputs(bfp->val, stderr);
165     fputs("'\n", stderr);
166     return(NULL);
167     }
168     zptr = (ZONE *)lep->data; /* add this subsurface */
169 greg 2.9 if (!idf_movobject(our_idf, param, zptr->plast))
170 greg 2.1 return(NULL);
171 greg 2.11 zptr->ntotal++;
172 greg 2.1 return(zptr);
173     }
174    
175 greg 2.9 /* Return field for vertices in the given object */
176 greg 2.3 static IDF_FIELD *
177 greg 2.9 get_vlist(IDF_OBJECT *param, const char *zname)
178 greg 2.2 {
179 greg 2.14 static const int itm_len = sizeof(IDF_FIELD)+6;
180     static char fld_buf[4*itm_len];
181     static char *next_fbp = fld_buf;
182     int i;
183     IDF_FIELD *res;
184 greg 2.11 /* check if subsurface */
185     if (!strcmp(param->pname, SUBSURF_PNAME)) {
186     if (zname != NULL) {
187     LUENT *lep = lu_find(&zonesurf_lut,
188     idf_getfield(param,SS_BASE_FLD)->val);
189     if (strcmp((*(ZONE *)lep->data).zname, zname))
190     return(NULL);
191     }
192 greg 2.13 res = idf_getfield(param, SS_VERT_FLD);
193     } else {
194     i = 0; /* check for surface type */
195     while (strcmp(surf_type[i].pname, param->pname))
196     if (surf_type[++i].pname == NULL)
197     return(NULL);
198    
199     if (zname != NULL) { /* matches specified zone? */
200     IDF_FIELD *fptr = idf_getfield(param, surf_type[i].zone_fld);
201     if (fptr == NULL || strcmp(fptr->val, zname))
202     return(NULL);
203     }
204     res = idf_getfield(param, surf_type[i].vert_fld);
205 greg 2.11 }
206 greg 2.13 if (!res->val[0]) { /* hack for missing #vert */
207     IDF_FIELD *fptr;
208     if (next_fbp >= fld_buf+sizeof(fld_buf))
209     next_fbp = fld_buf;
210     i = 0; /* count vertices */
211     for (fptr = res->next; fptr != NULL; fptr = fptr->next)
212     ++i;
213     if (i % 3)
214 greg 2.3 return(NULL);
215 greg 2.13 fptr = res->next;
216     res = (IDF_FIELD *)next_fbp; next_fbp += itm_len;
217     res->next = fptr;
218     res->rem = "";
219     sprintf(res->val, "%d", i/3);
220 greg 2.3 }
221 greg 2.13 return(res);
222 greg 2.2 }
223    
224 greg 2.11 /* Get/allocate surface polygon */
225     static SURFACE *
226     get_surface(IDF_FIELD *fptr)
227     {
228     SURFACE *surf;
229     int nv;
230     FVECT e1, e2, vc;
231     int i, j;
232     /* get number of vertices */
233     if (fptr == NULL || (nv = atoi(fptr->val)) < 3) {
234     fputs(progname, stderr);
235     fputs(": bad vertex count for surface!\n", stderr);
236     return(NULL);
237     }
238     surf = (SURFACE *)malloc(sizeof(SURFACE)+sizeof(FVECT)*(nv-3));
239     if (surf == NULL)
240     return(NULL);
241     surf->nv = nv;
242     for (i = nv; i--; ) /* reverse vertex order/orientation */
243     for (j = 0; j < 3; j++) {
244     fptr = fptr->next;
245     if (fptr == NULL) {
246     fputs(progname, stderr);
247     fputs(": missing vertex in get_surface()\n", stderr);
248     free(surf);
249     return(NULL);
250     }
251     if ((surf->vl[i][j] = atof(fptr->val)) == 0 &&
252     !isflt(fptr->val)) {
253     fputs(progname, stderr);
254     fputs(": bad vertex in get_surface()\n", stderr);
255     free(surf);
256     return(NULL);
257     }
258     }
259     /* compute area and normal */
260     surf->norm[0] = surf->norm[1] = surf->norm[2] = 0;
261     VSUB(e1, surf->vl[1], surf->vl[0]);
262     for (i = 2; i < nv; i++) {
263     VSUB(e2, surf->vl[i], surf->vl[0]);
264     fcross(vc, e1, e2);
265     surf->norm[0] += vc[0];
266     surf->norm[1] += vc[1];
267     surf->norm[2] += vc[2];
268     VCOPY(e1, e2);
269     }
270     surf->area = .5 * normalize(surf->norm);
271     if (surf->area == 0) {
272     fputs(progname, stderr);
273     fputs(": degenerate polygon in get_surface()\n", stderr);
274     free(surf);
275     return(NULL);
276     }
277     return(surf);
278     }
279    
280 greg 2.2 /* Convert surface to Radiance with modifier based on unique name */
281     static int
282 greg 2.9 rad_surface(IDF_OBJECT *param, FILE *ofp)
283 greg 2.2 {
284 greg 2.3 const char *sname = idf_getfield(param,NAME_FLD)->val;
285     IDF_FIELD *fptr = get_vlist(param, NULL);
286 greg 2.11 const char **fargs;
287     int nvert, i, j;
288 greg 2.2
289 greg 2.3 if (fptr == NULL || (nvert = atoi(fptr->val)) < 3) {
290 greg 2.2 fprintf(stderr, "%s: bad surface '%s'\n", progname, sname);
291     return(0);
292     }
293 greg 2.11 fargs = (const char **)malloc(sizeof(const char *)*3*nvert);
294     if (fargs == NULL)
295     return(0);
296 greg 2.2 fprintf(ofp, "\nvoid glow '%s'\n0\n0\n4 1 1 1 0\n", sname);
297     fprintf(ofp, "\n'%s' polygon 's_%s'\n0\n0\n%d\n", sname, sname, 3*nvert);
298 greg 2.11 for (j = nvert; j--; ) { /* get arguments in reverse */
299     for (i = 0; i < 3; i++) {
300 greg 2.2 fptr = fptr->next;
301 greg 2.3 if (fptr == NULL || !isflt(fptr->val)) {
302 greg 2.2 fprintf(stderr,
303 greg 2.3 "%s: missing/bad vertex for %s '%s'\n",
304     progname, param->pname, sname);
305 greg 2.2 return(0);
306     }
307 greg 2.11 fargs[3*j + i] = fptr->val;
308     }
309     }
310     for (j = 0; j < nvert; j++) { /* output reversed verts */
311     for (i = 0; i < 3; i++) {
312 greg 2.2 fputc('\t', ofp);
313 greg 2.11 fputs(fargs[3*j + i], ofp);
314 greg 2.2 }
315     fputc('\n', ofp);
316     }
317 greg 2.11 free(fargs);
318     return(!ferror(ofp));
319     }
320    
321     /* Convert subsurface to Radiance with modifier based on unique name */
322 greg 2.12 static double
323 greg 2.11 rad_subsurface(IDF_OBJECT *param, FILE *ofp)
324     {
325     const char *sname = idf_getfield(param,NAME_FLD)->val;
326 greg 2.14 SURFACE *surf = get_surface(get_vlist(param, NULL));
327 greg 2.12 double area;
328 greg 2.11 int i;
329    
330     if (surf == NULL) {
331     fprintf(stderr, "%s: bad subsurface '%s'\n", progname, sname);
332 greg 2.12 return(-1.);
333 greg 2.11 }
334     fprintf(ofp, "\nvoid glow '%s'\n0\n0\n4 1 1 1 0\n", sname);
335     fprintf(ofp, "\n'%s' polygon 'ss_%s'\n0\n0\n%d\n",
336     sname, sname, 3*surf->nv);
337     for (i = 0; i < surf->nv; i++) { /* offset vertices */
338     FVECT vert;
339     VSUM(vert, surf->vl[i], surf->norm, 2.*SURF_EPS);
340     fprintf(ofp, "\t%.12f %.12f %.12f\n", vert[0], vert[1], vert[2]);
341     }
342 greg 2.12 area = surf->area;
343 greg 2.11 free(surf);
344 greg 2.12 if (ferror(ofp))
345     return(-1.);
346     return(area);
347 greg 2.2 }
348    
349 greg 2.3 /* Start rcontrib process */
350 greg 2.2 static int
351     start_rcontrib(SUBPROC *pd, ZONE *zp)
352     {
353 greg 2.11 #define BASE_AC 6
354 greg 2.2 static char *base_av[BASE_AC] = {
355 greg 2.11 "rcontrib", "-V+", "-fff", "-h", "-x", "1"
356 greg 2.2 };
357     char cbuf[300];
358     char **av;
359     FILE *ofp;
360 greg 2.9 IDF_OBJECT *pptr;
361 greg 2.11 IDF_FIELD *fptr;
362 greg 2.12 int i, j, n;
363 greg 2.2 /* start oconv command */
364     sprintf(cbuf, "oconv - > '%s'", temp_octree);
365     if ((ofp = popen(cbuf, "w")) == NULL) {
366     fputs(progname, stderr);
367     fputs(": cannot open oconv process\n", stderr);
368     return(0);
369     }
370     /* allocate argument list */
371 greg 2.11 av = (char **)malloc(sizeof(char *)*(BASE_AC+4+2*(zp->ntotal)));
372 greg 2.2 if (av == NULL)
373     return(0);
374     for (i = 0; i < BASE_AC; i++)
375     av[i] = base_av[i];
376 greg 2.7 sprintf(cbuf, "%d", nsamps);
377 greg 2.2 av[i++] = "-c";
378 greg 2.11 av[i++] = cbuf; /* add modifiers & surfaces */
379     for (n = 0, pptr = zp->pfirst; n < zp->nsurf; n++, pptr = pptr->dnext) {
380     fptr = idf_getfield(pptr,NAME_FLD);
381 greg 2.3 if (fptr == NULL || !fptr->val[0]) {
382 greg 2.2 fputs(progname, stderr);
383 greg 2.9 fputs(": missing name for surface object\n", stderr);
384 greg 2.2 return(0);
385     }
386 greg 2.3 if (!rad_surface(pptr, ofp)) /* add surface to octree */
387     return(0);
388 greg 2.2 av[i++] = "-m";
389 greg 2.3 av[i++] = fptr->val;
390 greg 2.2 }
391 greg 2.11 /* now subsurfaces */
392 greg 2.12 if (zp->ntotal > zp->nsurf) {
393     if (zp->area_redu != NULL)
394     memset(zp->area_redu, 0, sizeof(float)*zp->nsurf);
395     else if ((zp->area_redu = (float *)calloc(zp->nsurf,
396     sizeof(float))) == NULL)
397     return(0);
398     }
399 greg 2.11 for ( ; n < zp->ntotal; n++, pptr = pptr->dnext) {
400 greg 2.12 double ssarea;
401     const char *bname;
402     IDF_OBJECT *pptr1;
403 greg 2.11 fptr = idf_getfield(pptr,NAME_FLD);
404     if (fptr == NULL || !fptr->val[0]) {
405     fputs(progname, stderr);
406     fputs(": missing name for subsurface object\n", stderr);
407     return(0);
408     }
409 greg 2.12 /* add subsurface to octree */
410     if ((ssarea = rad_subsurface(pptr, ofp)) < 0)
411 greg 2.11 return(0);
412 greg 2.12 /* mark area for subtraction */
413     bname = idf_getfield(pptr,SS_BASE_FLD)->val;
414     for (j = 0, pptr1 = zp->pfirst;
415     j < zp->nsurf; j++, pptr1 = pptr1->dnext)
416     if (!strcmp(idf_getfield(pptr1,NAME_FLD)->val, bname)) {
417     zp->area_redu[j] += ssarea;
418     break;
419     }
420 greg 2.11 av[i++] = "-m";
421     av[i++] = fptr->val;
422     }
423 greg 2.2 if (pclose(ofp) != 0) { /* finish oconv */
424     fputs(progname, stderr);
425     fputs(": error running oconv process\n", stderr);
426     return(0);
427     }
428     av[i++] = temp_octree; /* add final octree argument */
429     av[i] = NULL;
430 greg 2.8 if (open_process(pd, av) <= 0) { /* start process */
431 greg 2.2 fputs(progname, stderr);
432     fputs(": cannot start rcontrib process\n", stderr);
433     return(0);
434     }
435     free(av); /* all done -- clean up */
436     return(1);
437     #undef BASE_AC
438     }
439    
440 greg 2.3 /* Initialize polygon sampling */
441     static Vert2_list *
442 greg 2.11 init_poly(POLYSAMP *ps, IDF_FIELD *f0)
443 greg 2.3 {
444 greg 2.11 SURFACE *surf;
445     Vert2_list *vl2;
446     int i;
447     /* get 3-D polygon vertices */
448     if ((surf = get_surface(f0)) == NULL)
449     return(NULL);
450     vl2 = polyAlloc(surf->nv);
451 greg 2.3 if (vl2 == NULL)
452     return(NULL);
453 greg 2.11 /* create X & Y axes */
454     VCOPY(ps->sdir[2], surf->norm);
455     VSUB(ps->sdir[0], surf->vl[1], surf->vl[0]);
456     if (normalize(ps->sdir[0]) == 0)
457 greg 2.3 return(NULL);
458     fcross(ps->sdir[1], ps->sdir[2], ps->sdir[0]);
459     /* compute plane offset */
460 greg 2.11 ps->poff = DOT(surf->vl[0], ps->sdir[2]);
461 greg 2.3 /* assign 2-D vertices */
462 greg 2.11 for (i = 0; i < surf->nv; i++) {
463     vl2->v[i].mX = DOT(surf->vl[i], ps->sdir[0]);
464     vl2->v[i].mY = DOT(surf->vl[i], ps->sdir[1]);
465     }
466     ps->area_left = surf->area;
467     free(surf); /* it's ready! */
468     vl2->p = (void *)ps;
469 greg 2.3 return(vl2);
470     }
471    
472     /* Generate samples on 2-D triangle */
473     static int
474     sample_triangle(const Vert2_list *vl2, int a, int b, int c)
475     {
476     POLYSAMP *ps = (POLYSAMP *)vl2->p;
477     float *samp;
478     FVECT orig;
479     FVECT ab, ac;
480     double area;
481     int i, j, ns;
482     /* compute sampling axes */
483     for (i = 3; i--; ) {
484     orig[i] = vl2->v[a].mX*ps->sdir[0][i] +
485     vl2->v[a].mY*ps->sdir[1][i] +
486 greg 2.11 (ps->poff+SURF_EPS)*ps->sdir[2][i];
487 greg 2.3 ab[i] = (vl2->v[b].mX - vl2->v[a].mX)*ps->sdir[0][i] +
488     (vl2->v[b].mY - vl2->v[a].mY)*ps->sdir[1][i];
489     ac[i] = (vl2->v[c].mX - vl2->v[a].mX)*ps->sdir[0][i] +
490     (vl2->v[c].mY - vl2->v[a].mY)*ps->sdir[1][i];
491     }
492     /* compute number of samples to take */
493     area = .5*(vl2->v[a].mX*vl2->v[b].mY - vl2->v[b].mX*vl2->v[a].mY +
494     vl2->v[b].mX*vl2->v[c].mY - vl2->v[c].mX*vl2->v[b].mY +
495     vl2->v[c].mX*vl2->v[a].mY - vl2->v[a].mX*vl2->v[c].mY);
496     if (area < .0) {
497     fputs(progname, stderr);
498     fputs(": negative triangle area in sample_triangle()\n", stderr);
499     return(0);
500     }
501     if (area >= ps->area_left) {
502     ns = ps->samp_left;
503     ps->area_left = 0;
504     } else {
505     ns = (ps->samp_left*area/ps->area_left + .5);
506     ps->samp_left -= ns;
507     ps->area_left -= area;
508     }
509     if (ns <= 0) /* XXX should be error? */
510     return(1);
511     /* buffer sample rays */
512     samp = (float *)malloc(sizeof(float)*6*ns);
513     if (samp == NULL)
514     return(0);
515 greg 2.4 for (i = ns; i--; ) { /* stratified Monte Carlo sampling */
516 greg 2.3 double sv[4];
517 greg 2.5 FVECT dv;
518 greg 2.3 multisamp(sv, 4, (i+frandom())/(double)ns);
519     sv[0] *= sv[1] = sqrt(sv[1]);
520     sv[1] = 1. - sv[1];
521     for (j = 3; j--; )
522 greg 2.5 samp[i*6 + j] = orig[j] + sv[0]*ab[j] + sv[1]*ac[j];
523     sv[2] = sqrt(sv[2]);
524     sv[3] *= 2.*PI;
525     dv[0] = tcos(sv[3]) * sv[2];
526     dv[1] = tsin(sv[3]) * sv[2];
527     dv[2] = sqrt(1. - sv[2]*sv[2]);
528     for (j = 3; j--; )
529     samp[i*6 + 3 + j] = dv[0]*ps->sdir[0][j] +
530     dv[1]*ps->sdir[1][j] +
531     dv[2]*ps->sdir[2][j] ;
532 greg 2.3 }
533     /* send to our process */
534     writebuf(ps->wd, (char *)samp, sizeof(float)*6*ns);
535     free(samp); /* that's it! */
536     return(1);
537     }
538    
539     /* Sample the given surface */
540 greg 2.12 static double
541 greg 2.9 sample_surface(IDF_OBJECT *param, int wd)
542 greg 2.3 {
543     POLYSAMP psamp;
544 greg 2.12 double area;
545 greg 2.3 int nv;
546     Vert2_list *vlist2;
547     /* set up our polygon sampler */
548 greg 2.11 if ((vlist2 = init_poly(&psamp, get_vlist(param, NULL))) == NULL) {
549 greg 2.5 fprintf(stderr, "%s: bad polygon %s '%s'\n",
550 greg 2.3 progname, param->pname,
551     idf_getfield(param,NAME_FLD)->val);
552 greg 2.12 return(-1.);
553 greg 2.3 }
554 greg 2.7 psamp.samp_left = nsamps; /* assign samples & destination */
555 greg 2.3 psamp.wd = wd;
556 greg 2.11 /* hack for subsurface sampling */
557     psamp.poff += 2.*SURF_EPS * !strcmp(param->pname, SUBSURF_PNAME);
558 greg 2.12
559     area = psamp.area_left; /* remember starting surface area */
560 greg 2.3 /* sample each subtriangle */
561     if (!polyTriangulate(vlist2, &sample_triangle))
562 greg 2.12 return(-1.);
563 greg 2.3 polyFree(vlist2); /* clean up and return */
564 greg 2.12 return(area);
565 greg 2.3 }
566    
567 greg 2.2 /* Compute User View Factors using open rcontrib process */
568     static int
569 greg 2.3 compute_uvfs(SUBPROC *pd, ZONE *zp)
570 greg 2.2 {
571 greg 2.9 IDF_OBJECT *pptr, *pout, *pptr1;
572 greg 2.3 float *uvfa;
573     char uvfbuf[24];
574     int n, m;
575 greg 2.9 /* create output object */
576     pout = idf_newobject(our_idf, UVF_PNAME,
577 greg 2.11 " ! computed by Radiance\n ", our_idf->plast);
578 greg 2.3 if (pout == NULL) {
579     fputs(progname, stderr);
580 greg 2.9 fputs(": cannot create new IDF object\n", stderr);
581 greg 2.3 return(0);
582     }
583     if (!idf_addfield(pout, zp->zname,
584     " !- Zone Name\n ")) {
585     fputs(progname, stderr);
586     fputs(": cannot add zone name field\n", stderr);
587     return(0);
588     }
589     /* allocate read buffer */
590 greg 2.11 uvfa = (float *)malloc(sizeof(float)*3*zp->ntotal);
591 greg 2.3 if (uvfa == NULL)
592     return(0);
593     /* UVFs from each surface */
594 greg 2.11 for (n = 0, pptr = zp->pfirst; n < zp->ntotal; n++, pptr = pptr->dnext) {
595 greg 2.4 double vfsum = 0;
596 greg 2.12 double adj_factor;
597 greg 2.3 /* send samples to rcontrib */
598 greg 2.12 if ((adj_factor = sample_surface(pptr, pd->w)) < 0)
599 greg 2.3 return(0);
600 greg 2.12 if (zp->area_redu == NULL)
601     adj_factor = 1.;
602     else /* comp. for subsurface area */
603     adj_factor /= adj_factor - zp->area_redu[n];
604 greg 2.3 /* read results */
605 greg 2.11 if (readbuf(pd->r, (char *)uvfa, sizeof(float)*3*zp->ntotal) !=
606     sizeof(float)*3*zp->ntotal) {
607 greg 2.3 fputs(progname, stderr);
608     fputs(": read error from rcontrib process\n", stderr);
609     return(0);
610     }
611     /* append UVF fields */
612     for (m = 0, pptr1 = zp->pfirst;
613 greg 2.11 m < zp->ntotal; m++, pptr1 = pptr1->dnext) {
614 greg 2.12 const double uvf = uvfa[3*m + 1] * adj_factor;
615     vfsum += uvf;
616     if (pptr1 == pptr && uvf > .001)
617 greg 2.11 fprintf(stderr,
618 greg 2.5 "%s: warning - non-zero self-VF (%.1f%%) for surface '%s'\n",
619 greg 2.12 progname, 100.*uvf,
620 greg 2.3 idf_getfield(pptr,NAME_FLD)->val);
621 greg 2.12 sprintf(uvfbuf, "%.4f", uvf);
622 greg 2.3 if (!idf_addfield(pout,
623     idf_getfield(pptr,NAME_FLD)->val, NULL) ||
624     !idf_addfield(pout,
625     idf_getfield(pptr1,NAME_FLD)->val, NULL) ||
626     !idf_addfield(pout, uvfbuf,
627 greg 2.11 (n+m < 2*zp->ntotal-2) ?
628     "\n " : "\n\n")) {
629 greg 2.3 fputs(progname, stderr);
630     fputs(": error adding UVF fields\n", stderr);
631     return(0);
632     }
633     }
634 greg 2.4 if (vfsum < 0.95)
635     fprintf(stderr,
636     "%s: warning - missing %.1f%% of energy from surface '%s'\n",
637     progname, 100.*(1.-vfsum),
638     idf_getfield(pptr,NAME_FLD)->val);
639 greg 2.3 }
640     free(uvfa); /* clean up and return */
641     return(1);
642 greg 2.2 }
643    
644 greg 2.1 /* Compute zone User View Factors */
645     static int
646     compute_zones(void)
647     {
648 greg 2.2 ZONE *zptr;
649     /* temporary octree name */
650 greg 2.5 mktemp(strcpy(temp_octree, TEMPLATE));
651 greg 2.2 /* compute each zone */
652     for (zptr = zone_list; zptr != NULL; zptr = zptr->next) {
653     SUBPROC rcproc;
654     /* start rcontrib process */
655     if (!start_rcontrib(&rcproc, zptr))
656     return(0);
657     /* compute+add view factors */
658     if (!compute_uvfs(&rcproc, zptr))
659     return(0);
660     if (close_process(&rcproc) != 0) {
661     fputs(progname, stderr);
662     fputs(": bad return status from rcontrib\n", stderr);
663     return(0);
664     }
665     }
666     unlink(temp_octree); /* remove octree file */
667     return(1);
668 greg 2.1 }
669    
670     /* Load IDF and compute User View Factors */
671     int
672     main(int argc, char *argv[])
673     {
674 greg 2.2 int incl_comments = 1;
675 greg 2.1 char *origIDF, *revIDF;
676 greg 2.9 IDF_OBJECT *pptr;
677 greg 2.1 int i;
678    
679 greg 2.7 progname = *argv++; argc--; /* get options if any */
680     while (argc > 1 && argv[0][0] == '-')
681     switch (argv[0][1]) {
682     case 'c': /* elide comments */
683     incl_comments = -1; /* header only */
684     argv++; argc--;
685     continue;
686     case 's': /* samples */
687     nsamps = 1000*atoi(*++argv);
688     argv++; argc -= 2;
689     continue;
690     default:
691     fputs(progname, stderr);
692     fputs(": unknown option '", stderr);
693     fputs(argv[0], stderr);
694     fputs("'\n", stderr);
695     goto userr;
696     }
697     if ((argc < 1) | (argc > 2))
698     goto userr;
699     origIDF = argv[0];
700     revIDF = (argc == 1) ? argv[0] : argv[1];
701 greg 2.1 /* load Input Data File */
702     our_idf = idf_load(origIDF);
703     if (our_idf == NULL) {
704 greg 2.2 fputs(progname, stderr);
705 greg 2.1 fputs(": cannot load IDF '", stderr);
706     fputs(origIDF, stderr);
707     fputs("'\n", stderr);
708     return(1);
709     }
710     /* remove existing UVFs */
711 greg 2.9 if ((pptr = idf_getobject(our_idf, UVF_PNAME)) != NULL) {
712     IDF_OBJECT *pnext;
713 greg 2.2 fputs(progname, stderr);
714 greg 2.1 fputs(": removing previous User View Factors\n", stderr);
715     do {
716     pnext = pptr->pnext;
717 greg 2.9 idf_delobject(our_idf, pptr);
718 greg 2.1 } while (pnext != NULL);
719     }
720     /* add to header */
721     if (our_idf->hrem == NULL ||
722     (i = strlen(our_idf->hrem)-strlen(ADD_HEADER)) < 0 ||
723     strcmp(our_idf->hrem+i, ADD_HEADER))
724     idf_add2hdr(our_idf, ADD_HEADER);
725     /* gather zone surfaces */
726     for (i = 0; surf_type[i].pname != NULL; i++)
727 greg 2.9 for (pptr = idf_getobject(our_idf, surf_type[i].pname);
728 greg 2.1 pptr != NULL; pptr = pptr->pnext) {
729     IDF_FIELD *fptr = idf_getfield(pptr,
730     surf_type[i].zone_fld);
731     if (fptr == NULL) {
732 greg 2.2 fputs(progname, stderr);
733     fputs(": warning - missing zone field\n", stderr);
734 greg 2.1 continue;
735     }
736 greg 2.3 if (add2zone(pptr, fptr->val) == NULL)
737 greg 2.1 return(1);
738     }
739 greg 2.11 /* add subsurfaces */
740     for (pptr = idf_getobject(our_idf, SUBSURF_PNAME);
741     pptr != NULL; pptr = pptr->pnext)
742     if (add_subsurf(pptr) == NULL)
743     return(1);
744 greg 2.1 /* run rcontrib on each zone */
745     if (!compute_zones())
746     return(1);
747     /* write out modified IDF */
748 greg 2.2 if (!idf_write(our_idf, revIDF, incl_comments)) {
749     fputs(progname, stderr);
750 greg 2.1 fputs(": error writing IDF '", stderr);
751     fputs(revIDF, stderr);
752     fputs("'\n", stderr);
753     return(1);
754     }
755     return(0); /* finito! */
756 greg 2.7 userr:
757     fputs("Usage: ", stderr);
758     fputs(progname, stderr);
759     fputs(" [-c][-s Ksamps] Model.idf [Revised.idf]\n", stderr);
760     return(1);
761 greg 2.1 }