ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/util/eplus_adduvf.c
Revision: 2.18
Committed: Sun Mar 6 01:13:18 2016 UTC (8 years, 1 month ago) by schorsch
Content type: text/plain
Branch: MAIN
CVS Tags: rad5R2, rad5R1
Changes since 2.17: +2 -1 lines
Log Message:
Prepare for SCons build on Win32 and Win64

File Contents

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