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

Comparing ray/src/util/eplus_adduvf.c (file contents):
Revision 2.10 by greg, Thu Feb 13 17:33:37 2014 UTC vs.
Revision 2.12 by greg, Fri Feb 21 13:17:45 2014 UTC

# Line 19 | Line 19 | static const char RCSid[] = "$Id$";
19   #define NSAMPLES        80000                   /* default number of samples */
20   #endif
21  
22 + #define SURF_EPS        0.0005                  /* surface testing epsilon */
23 +
24   char            *progname;                      /* global argv[0] */
25  
26   int             nsamps = NSAMPLES;              /* number of samples to use */
# Line 28 | Line 30 | char           temp_octree[128];               /* temporary octree */
30   const char      UVF_PNAME[] =
31                          "ZoneProperty:UserViewFactors:bySurfaceName";
32  
33 + const char      SUBSURF_PNAME[] =
34 +                        "FenestrationSurface:Detailed";
35 +
36   const char      ADD_HEADER[] =
37                          "\n!+++ User View Factors computed by Radiance +++!\n\n";
38  
39   #define NAME_FLD        1                       /* name field always first? */
40  
41 + #define SS_BASE_FLD     4                       /* subsurface base surface */
42 + #define SS_VERT_FLD     10                      /* subsurface vertex count */
43 +
44   typedef struct {
45          const char      *pname;                 /* object type name */
46          short           zone_fld;               /* zone field index */
# Line 51 | Line 59 | typedef struct s_zone {
59          const char      *zname;                 /* zone name */
60          struct s_zone   *next;                  /* next zone in list */
61          int             nsurf;                  /* surface count */
62 +        int             ntotal;                 /* surfaces+subsurfaces */
63          IDF_OBJECT      *pfirst;                /* first matching object */
64 <        IDF_OBJECT      *plast;                 /* last matching object */
64 >        IDF_OBJECT      *plast;                 /* last before subsurfaces */
65 >        float           *area_redu;             /* subsurface area per surf. */
66   } ZONE;                 /* a list of collected zone surfaces */
67  
68   ZONE            *zone_list = NULL;      /* our list of zones */
69  
70 + LUTAB           zonesurf_lut =          /* zone surface lookup table */
71 +                        LU_SINIT(NULL,NULL);
72 +
73   IDF_LOADED      *our_idf = NULL;        /* loaded/modified IDF */
74  
75   typedef struct {
76 +        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          FVECT           sdir[3];        /* UVW unit sampling vectors */
84          double          poff;           /* W-offset for plane of polygon */
85          double          area_left;      /* area left to sample */
# Line 78 | Line 98 | new_zone(const char *zname, IDF_OBJECT *param)
98          znew->zname = zname;                    /* assumes static */
99          znew->next = zone_list;
100          znew->pfirst = znew->plast = param;
101 <        znew->nsurf = 1;
101 >        znew->ntotal = znew->nsurf = 1;
102 >        znew->area_redu = NULL;
103          return(zone_list = znew);
104   }
105  
# Line 86 | Line 107 | new_zone(const char *zname, IDF_OBJECT *param)
107   static ZONE *
108   add2zone(IDF_OBJECT *param, const char *zname)
109   {
110 <        ZONE    *zptr;
110 >        IDF_FIELD       *nfp = idf_getfield(param, NAME_FLD);
111 >        ZONE            *zptr;
112 >        LUENT           *lep;
113  
114 +        if (nfp == NULL) {
115 +                fputs(progname, stderr);
116 +                fputs(": surface missing name field!\n", stderr);
117 +                return(NULL);
118 +        }
119          for (zptr = zone_list; zptr != NULL; zptr = zptr->next)
120                  if (!strcmp(zptr->zname, zname))
121                          break;
122 <        if (zptr == NULL)
123 <                return(new_zone(zname, param));
124 <                                                /* keep surfaces together */
122 >        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          if (!idf_movobject(our_idf, param, zptr->plast))
170                  return(NULL);
171 <        zptr->plast = param;
100 <        zptr->nsurf++;
171 >        zptr->ntotal++;
172          return(zptr);
173   }
174  
# Line 106 | Line 177 | static IDF_FIELD *
177   get_vlist(IDF_OBJECT *param, const char *zname)
178   {
179          int             i = 0;
180 <        IDF_FIELD       *fptr;
180 >                                                /* check if subsurface */
181 >        if (!strcmp(param->pname, SUBSURF_PNAME)) {
182 >                if (zname != NULL) {
183 >                        LUENT   *lep = lu_find(&zonesurf_lut,
184 >                                        idf_getfield(param,SS_BASE_FLD)->val);
185 >                        if (strcmp((*(ZONE *)lep->data).zname, zname))
186 >                                return(NULL);
187 >                }
188 >                return(idf_getfield(param, SS_VERT_FLD));
189 >        }
190                                                  /* check for surface type */
191          while (strcmp(surf_type[i].pname, param->pname))
192                  if (surf_type[++i].pname == NULL)
193                          return(NULL);
194  
195          if (zname != NULL) {                    /* matches specified zone? */
196 <                fptr = idf_getfield(param, surf_type[i].zone_fld);
196 >                IDF_FIELD       *fptr = idf_getfield(param, surf_type[i].zone_fld);
197                  if (fptr == NULL || strcmp(fptr->val, zname))
198                          return(NULL);
199          }
# Line 121 | Line 201 | get_vlist(IDF_OBJECT *param, const char *zname)
201          return(idf_getfield(param, surf_type[i].vert_fld));
202   }
203  
204 + /* Get/allocate surface polygon */
205 + static SURFACE *
206 + get_surface(IDF_FIELD *fptr)
207 + {
208 +        SURFACE *surf;
209 +        int     nv;
210 +        FVECT   e1, e2, vc;
211 +        int     i, j;
212 +                                        /* get number of vertices */
213 +        if (fptr == NULL || (nv = atoi(fptr->val)) < 3) {
214 +                fputs(progname, stderr);
215 +                fputs(": bad vertex count for surface!\n", stderr);
216 +                return(NULL);
217 +        }
218 +        surf = (SURFACE *)malloc(sizeof(SURFACE)+sizeof(FVECT)*(nv-3));
219 +        if (surf == NULL)
220 +                return(NULL);
221 +        surf->nv = nv;
222 +        for (i = nv; i--; )             /* reverse vertex order/orientation */
223 +                for (j = 0; j < 3; j++) {
224 +                        fptr = fptr->next;
225 +                        if (fptr == NULL) {
226 +                                fputs(progname, stderr);
227 +                                fputs(": missing vertex in get_surface()\n", stderr);
228 +                                free(surf);
229 +                                return(NULL);
230 +                        }
231 +                        if ((surf->vl[i][j] = atof(fptr->val)) == 0 &&
232 +                                                        !isflt(fptr->val)) {
233 +                                fputs(progname, stderr);
234 +                                fputs(": bad vertex in get_surface()\n", stderr);
235 +                                free(surf);
236 +                                return(NULL);
237 +                        }
238 +                }
239 +                                        /* compute area and normal */
240 +        surf->norm[0] = surf->norm[1] = surf->norm[2] = 0;
241 +        VSUB(e1, surf->vl[1], surf->vl[0]);
242 +        for (i = 2; i < nv; i++) {
243 +                VSUB(e2, surf->vl[i], surf->vl[0]);
244 +                fcross(vc, e1, e2);
245 +                surf->norm[0] += vc[0];
246 +                surf->norm[1] += vc[1];
247 +                surf->norm[2] += vc[2];
248 +                VCOPY(e1, e2);
249 +        }
250 +        surf->area = .5 * normalize(surf->norm);
251 +        if (surf->area == 0) {
252 +                fputs(progname, stderr);
253 +                fputs(": degenerate polygon in get_surface()\n", stderr);
254 +                free(surf);
255 +                return(NULL);
256 +        }
257 +        return(surf);
258 + }
259 +
260   /* Convert surface to Radiance with modifier based on unique name */
261   static int
262   rad_surface(IDF_OBJECT *param, FILE *ofp)
263   {
264          const char      *sname = idf_getfield(param,NAME_FLD)->val;
265          IDF_FIELD       *fptr = get_vlist(param, NULL);
266 <        int             nvert, i;
266 >        const char      **fargs;
267 >        int             nvert, i, j;
268  
269          if (fptr == NULL || (nvert = atoi(fptr->val)) < 3) {
270                  fprintf(stderr, "%s: bad surface '%s'\n", progname, sname);
271                  return(0);
272          }
273 +        fargs = (const char **)malloc(sizeof(const char *)*3*nvert);
274 +        if (fargs == NULL)
275 +                return(0);
276          fprintf(ofp, "\nvoid glow '%s'\n0\n0\n4 1 1 1 0\n", sname);
277          fprintf(ofp, "\n'%s' polygon 's_%s'\n0\n0\n%d\n", sname, sname, 3*nvert);
278 <        while (nvert--) {
279 <                for (i = 3; i--; ) {
278 >        for (j = nvert; j--; ) {                /* get arguments in reverse */
279 >                for (i = 0; i < 3; i++) {
280                          fptr = fptr->next;
281                          if (fptr == NULL || !isflt(fptr->val)) {
282                                  fprintf(stderr,
# Line 144 | Line 284 | rad_surface(IDF_OBJECT *param, FILE *ofp)
284                                                  progname, param->pname, sname);
285                                  return(0);
286                          }
287 +                        fargs[3*j + i] = fptr->val;
288 +                }
289 +        }
290 +        for (j = 0; j < nvert; j++) {           /* output reversed verts */
291 +                for (i = 0; i < 3; i++) {
292                          fputc('\t', ofp);
293 <                        fputs(fptr->val, ofp);
293 >                        fputs(fargs[3*j + i], ofp);
294                  }
295                  fputc('\n', ofp);
296          }
297 +        free(fargs);
298          return(!ferror(ofp));
299   }
300  
301 + /* Convert subsurface to Radiance with modifier based on unique name */
302 + static double
303 + rad_subsurface(IDF_OBJECT *param, FILE *ofp)
304 + {
305 +        const char      *sname = idf_getfield(param,NAME_FLD)->val;
306 +        SURFACE         *surf = get_surface(idf_getfield(param,SS_VERT_FLD));
307 +        double          area;
308 +        int             i;
309 +
310 +        if (surf == NULL) {
311 +                fprintf(stderr, "%s: bad subsurface '%s'\n", progname, sname);
312 +                return(-1.);
313 +        }
314 +        fprintf(ofp, "\nvoid glow '%s'\n0\n0\n4 1 1 1 0\n", sname);
315 +        fprintf(ofp, "\n'%s' polygon 'ss_%s'\n0\n0\n%d\n",
316 +                                                sname, sname, 3*surf->nv);
317 +        for (i = 0; i < surf->nv; i++) {        /* offset vertices */
318 +                FVECT   vert;
319 +                VSUM(vert, surf->vl[i], surf->norm, 2.*SURF_EPS);
320 +                fprintf(ofp, "\t%.12f %.12f %.12f\n", vert[0], vert[1], vert[2]);
321 +        }
322 +        area = surf->area;
323 +        free(surf);
324 +        if (ferror(ofp))
325 +                return(-1.);
326 +        return(area);
327 + }
328 +
329   /* Start rcontrib process */
330   static int
331   start_rcontrib(SUBPROC *pd, ZONE *zp)
332   {
333 < #define BASE_AC         5
333 > #define BASE_AC         6
334          static char     *base_av[BASE_AC] = {
335 <                                "rcontrib", "-fff", "-h", "-x", "1"
335 >                                "rcontrib", "-V+", "-fff", "-h", "-x", "1"
336                          };
337          char            cbuf[300];
338          char            **av;
339          FILE            *ofp;
340          IDF_OBJECT      *pptr;
341 <        int             i, n;
341 >        IDF_FIELD       *fptr;
342 >        int             i, j, n;
343                                                  /* start oconv command */
344          sprintf(cbuf, "oconv - > '%s'", temp_octree);
345          if ((ofp = popen(cbuf, "w")) == NULL) {
# Line 173 | Line 348 | start_rcontrib(SUBPROC *pd, ZONE *zp)
348                  return(0);
349          }
350                                                  /* allocate argument list */
351 <        av = (char **)malloc(sizeof(char *)*(BASE_AC+4+2*zp->nsurf));
351 >        av = (char **)malloc(sizeof(char *)*(BASE_AC+4+2*(zp->ntotal)));
352          if (av == NULL)
353                  return(0);
354          for (i = 0; i < BASE_AC; i++)
355                  av[i] = base_av[i];
356          sprintf(cbuf, "%d", nsamps);
357          av[i++] = "-c";
358 <        av[i++] = cbuf;                         /* add modifier arguments */
359 <        for (n = zp->nsurf, pptr = zp->pfirst; n--; pptr = pptr->dnext) {
360 <                IDF_FIELD       *fptr = idf_getfield(pptr,NAME_FLD);
358 >        av[i++] = cbuf;                         /* add modifiers & surfaces */
359 >        for (n = 0, pptr = zp->pfirst; n < zp->nsurf; n++, pptr = pptr->dnext) {
360 >                fptr = idf_getfield(pptr,NAME_FLD);
361                  if (fptr == NULL || !fptr->val[0]) {
362                          fputs(progname, stderr);
363                          fputs(": missing name for surface object\n", stderr);
# Line 193 | Line 368 | start_rcontrib(SUBPROC *pd, ZONE *zp)
368                  av[i++] = "-m";
369                  av[i++] = fptr->val;
370          }
371 +                                                /* now subsurfaces */
372 +        if (zp->ntotal > zp->nsurf) {
373 +                if (zp->area_redu != NULL)
374 +                        memset(zp->area_redu, 0, sizeof(float)*zp->nsurf);
375 +                else if ((zp->area_redu = (float *)calloc(zp->nsurf,
376 +                                                sizeof(float))) == NULL)
377 +                        return(0);
378 +        }
379 +        for ( ; n < zp->ntotal; n++, pptr = pptr->dnext) {
380 +                double          ssarea;
381 +                const char      *bname;
382 +                IDF_OBJECT      *pptr1;
383 +                fptr = idf_getfield(pptr,NAME_FLD);
384 +                if (fptr == NULL || !fptr->val[0]) {
385 +                        fputs(progname, stderr);
386 +                        fputs(": missing name for subsurface object\n", stderr);
387 +                        return(0);
388 +                }
389 +                                                /* add subsurface to octree */
390 +                if ((ssarea = rad_subsurface(pptr, ofp)) < 0)
391 +                        return(0);
392 +                                                /* mark area for subtraction */
393 +                bname = idf_getfield(pptr,SS_BASE_FLD)->val;
394 +                for (j = 0, pptr1 = zp->pfirst;
395 +                                j < zp->nsurf; j++, pptr1 = pptr1->dnext)
396 +                        if (!strcmp(idf_getfield(pptr1,NAME_FLD)->val, bname)) {
397 +                                zp->area_redu[j] += ssarea;
398 +                                break;
399 +                        }
400 +                av[i++] = "-m";
401 +                av[i++] = fptr->val;
402 +        }
403          if (pclose(ofp) != 0) {                 /* finish oconv */
404                  fputs(progname, stderr);
405                  fputs(": error running oconv process\n", stderr);
# Line 212 | Line 419 | start_rcontrib(SUBPROC *pd, ZONE *zp)
419  
420   /* Initialize polygon sampling */
421   static Vert2_list *
422 < init_poly(POLYSAMP *ps, IDF_FIELD *f0, int nv)
422 > init_poly(POLYSAMP *ps, IDF_FIELD *f0)
423   {
424 <        IDF_FIELD       *fptr = f0;
425 <        int             i, j;
426 <        FVECT           *vl3, e1, e2, vc;
427 <        Vert2_list      *vl2 = polyAlloc(nv);
428 <
424 >        SURFACE         *surf;
425 >        Vert2_list      *vl2;
426 >        int             i;
427 >                                        /* get 3-D polygon vertices */
428 >        if ((surf = get_surface(f0)) == NULL)
429 >                return(NULL);
430 >        vl2 = polyAlloc(surf->nv);
431          if (vl2 == NULL)
432                  return(NULL);
224 –        vl2->p = ps;
225 –                                        /* get 3-D vertices */
226 –        vl3 = (FVECT *)malloc(sizeof(FVECT)*nv);
227 –        if (vl3 == NULL)
228 –                return(NULL);
229 –        for (i = nv; i--; )             /* reverse vertex ordering */
230 –                for (j = 0; j < 3; j++) {
231 –                        if (fptr == NULL) {
232 –                                fputs(progname, stderr);
233 –                                fputs(": missing vertex in init_poly()\n", stderr);
234 –                                return(NULL);
235 –                        }
236 –                        vl3[i][j] = atof(fptr->val);
237 –                        fptr = fptr->next;
238 –                }
239 –                                        /* compute area and normal */
240 –        ps->sdir[2][0] = ps->sdir[2][1] = ps->sdir[2][2] = 0;
241 –        VSUB(e1, vl3[1], vl3[0]);
242 –        for (i = 2; i < nv; i++) {
243 –                VSUB(e2, vl3[i], vl3[0]);
244 –                fcross(vc, e1, e2);
245 –                ps->sdir[2][0] += vc[0];
246 –                ps->sdir[2][1] += vc[1];
247 –                ps->sdir[2][2] += vc[2];
248 –                VCOPY(e1, e2);
249 –        }
250 –        ps->area_left = .5 * normalize(ps->sdir[2]);
251 –        if (ps->area_left == .0) {
252 –                fputs(progname, stderr);
253 –                fputs(": degenerate polygon in init_poly()\n", stderr);
254 –                return(0);
255 –        }
433                                          /* create X & Y axes */
434 <        VCOPY(ps->sdir[0], e1);
435 <        normalize(ps->sdir[0]);
434 >        VCOPY(ps->sdir[2], surf->norm);
435 >        VSUB(ps->sdir[0], surf->vl[1], surf->vl[0]);
436 >        if (normalize(ps->sdir[0]) == 0)
437 >                return(NULL);
438          fcross(ps->sdir[1], ps->sdir[2], ps->sdir[0]);
439                                          /* compute plane offset */
440 <        ps->poff = DOT(vl3[0], ps->sdir[2]);
440 >        ps->poff = DOT(surf->vl[0], ps->sdir[2]);
441                                          /* assign 2-D vertices */
442 <        for (i = 0; i < nv; i++) {
443 <                vl2->v[i].mX = DOT(vl3[i], ps->sdir[0]);
444 <                vl2->v[i].mY = DOT(vl3[i], ps->sdir[1]);
442 >        for (i = 0; i < surf->nv; i++) {
443 >                vl2->v[i].mX = DOT(surf->vl[i], ps->sdir[0]);
444 >                vl2->v[i].mY = DOT(surf->vl[i], ps->sdir[1]);
445          }
446 <        free(vl3);                      /* it's ready! */
446 >        ps->area_left = surf->area;
447 >        free(surf);                     /* it's ready! */
448 >        vl2->p = (void *)ps;
449          return(vl2);
450   }
451  
# Line 282 | Line 463 | sample_triangle(const Vert2_list *vl2, int a, int b, i
463          for (i = 3; i--; ) {
464                  orig[i] = vl2->v[a].mX*ps->sdir[0][i] +
465                                  vl2->v[a].mY*ps->sdir[1][i] +
466 <                                (ps->poff+.001)*ps->sdir[2][i];
466 >                                (ps->poff+SURF_EPS)*ps->sdir[2][i];
467                  ab[i] = (vl2->v[b].mX - vl2->v[a].mX)*ps->sdir[0][i] +
468                                  (vl2->v[b].mY - vl2->v[a].mY)*ps->sdir[1][i];
469                  ac[i] = (vl2->v[c].mX - vl2->v[a].mX)*ps->sdir[0][i] +
# Line 336 | Line 517 | sample_triangle(const Vert2_list *vl2, int a, int b, i
517   }
518  
519   /* Sample the given surface */
520 < static int
520 > static double
521   sample_surface(IDF_OBJECT *param, int wd)
522   {
342 –        IDF_FIELD       *fptr = get_vlist(param, NULL);
523          POLYSAMP        psamp;
524 +        double          area;
525          int             nv;
526          Vert2_list      *vlist2;
527                                          /* set up our polygon sampler */
528 <        if (fptr == NULL || (nv = atoi(fptr->val)) < 3 ||
348 <                        (vlist2 = init_poly(&psamp, fptr->next, nv)) == NULL) {
528 >        if ((vlist2 = init_poly(&psamp, get_vlist(param, NULL))) == NULL) {
529                  fprintf(stderr, "%s: bad polygon %s '%s'\n",
530                                  progname, param->pname,
531                                  idf_getfield(param,NAME_FLD)->val);
532 <                return(0);
532 >                return(-1.);
533          }
534          psamp.samp_left = nsamps;       /* assign samples & destination */
535          psamp.wd = wd;
536 +                                        /* hack for subsurface sampling */
537 +        psamp.poff += 2.*SURF_EPS * !strcmp(param->pname, SUBSURF_PNAME);
538 +
539 +        area = psamp.area_left;         /* remember starting surface area */
540                                          /* sample each subtriangle */
541          if (!polyTriangulate(vlist2, &sample_triangle))
542 <                return(0);
542 >                return(-1.);
543          polyFree(vlist2);               /* clean up and return */
544 <        return(1);
544 >        return(area);
545   }
546  
547   /* Compute User View Factors using open rcontrib process */
# Line 370 | Line 554 | compute_uvfs(SUBPROC *pd, ZONE *zp)
554          int             n, m;
555                                                  /* create output object */
556          pout = idf_newobject(our_idf, UVF_PNAME,
557 <                        "    ! computed by Radiance\n        ", zp->plast);
557 >                        "    ! computed by Radiance\n        ", our_idf->plast);
558          if (pout == NULL) {
559                  fputs(progname, stderr);
560                  fputs(": cannot create new IDF object\n", stderr);
# Line 383 | Line 567 | compute_uvfs(SUBPROC *pd, ZONE *zp)
567                  return(0);
568          }
569                                                  /* allocate read buffer */
570 <        uvfa = (float *)malloc(sizeof(float)*3*zp->nsurf);
570 >        uvfa = (float *)malloc(sizeof(float)*3*zp->ntotal);
571          if (uvfa == NULL)
572                  return(0);
573                                                  /* UVFs from each surface */
574 <        for (n = zp->nsurf, pptr = zp->pfirst; n--; pptr = pptr->dnext) {
574 >        for (n = 0, pptr = zp->pfirst; n < zp->ntotal; n++, pptr = pptr->dnext) {
575                  double  vfsum = 0;
576 +                double  adj_factor;
577                                                  /* send samples to rcontrib */
578 <                if (!sample_surface(pptr, pd->w))
578 >                if ((adj_factor = sample_surface(pptr, pd->w)) < 0)
579                          return(0);
580 +                if (zp->area_redu == NULL)
581 +                        adj_factor = 1.;
582 +                else                            /* comp. for subsurface area */
583 +                        adj_factor /= adj_factor - zp->area_redu[n];
584                                                  /* read results */
585 <                if (readbuf(pd->r, (char *)uvfa, sizeof(float)*3*zp->nsurf) !=
586 <                                sizeof(float)*3*zp->nsurf) {
585 >                if (readbuf(pd->r, (char *)uvfa, sizeof(float)*3*zp->ntotal) !=
586 >                                sizeof(float)*3*zp->ntotal) {
587                          fputs(progname, stderr);
588                          fputs(": read error from rcontrib process\n", stderr);
589                          return(0);
590                  }
591                                                  /* append UVF fields */
592                  for (m = 0, pptr1 = zp->pfirst;
593 <                                m < zp->nsurf; m++, pptr1 = pptr1->dnext) {
594 <                        vfsum += uvfa[3*m + 1];
595 <                        if (pptr1 == pptr) {
596 <                                if (uvfa[3*m + 1] > .001)
597 <                                        fprintf(stderr,
593 >                                m < zp->ntotal; m++, pptr1 = pptr1->dnext) {
594 >                        const double    uvf = uvfa[3*m + 1] * adj_factor;
595 >                        vfsum += uvf;
596 >                        if (pptr1 == pptr && uvf > .001)
597 >                                fprintf(stderr,
598                  "%s: warning - non-zero self-VF (%.1f%%) for surface '%s'\n",
599 <                                                progname, 100.*uvfa[3*m + 1],
599 >                                                progname, 100.*uvf,
600                                                  idf_getfield(pptr,NAME_FLD)->val);
601 <                                continue;       /* don't record self-factor */
413 <                        }
414 <                        sprintf(uvfbuf, "%.4f", uvfa[3*m + 1]);
601 >                        sprintf(uvfbuf, "%.4f", uvf);
602                          if (!idf_addfield(pout,
603                                          idf_getfield(pptr,NAME_FLD)->val, NULL) ||
604                                  !idf_addfield(pout,
605                                          idf_getfield(pptr1,NAME_FLD)->val, NULL) ||
606                                  !idf_addfield(pout, uvfbuf,
607 <                                                (n || m < zp->nsurf-2) ?
608 <                                                        "\n        " : "\n\n")) {
607 >                                                (n+m < 2*zp->ntotal-2) ?
608 >                                                "\n        " : "\n\n")) {
609                                  fputs(progname, stderr);
610                                  fputs(": error adding UVF fields\n", stderr);
611                                  return(0);
# Line 529 | Line 716 | main(int argc, char *argv[])
716                          if (add2zone(pptr, fptr->val) == NULL)
717                                  return(1);
718                  }
719 +                                                /* add subsurfaces */
720 +        for (pptr = idf_getobject(our_idf, SUBSURF_PNAME);
721 +                        pptr != NULL; pptr = pptr->pnext)
722 +                if (add_subsurf(pptr) == NULL)
723 +                        return(1);
724                                                  /* run rcontrib on each zone */
725          if (!compute_zones())
726                  return(1);

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines