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

Comparing ray/src/util/glarendx.c (file contents):
Revision 2.2 by greg, Tue Apr 28 09:11:00 1992 UTC vs.
Revision 2.9 by greg, Wed Nov 19 16:21:28 2003 UTC

# Line 1 | Line 1
1 – /* Copyright (c) 1991 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   * Compute Glare Index given by program name or -t option:
6   *
7 + *      dgi -           Daylight Glare Index
8 + *      brs_gi -        Building Research Station Glare Index (Petherbridge
9 + *                                                             & Hopkinson)
10 + *      ugr -           Unified Glare Rating System (Fischer)
11   *      guth_dgr -      Guth discomfort glare rating
12   *      guth_vcp -      Guth visual comfort probability
13   *      cie_cgi -       CIE Glare Index (1983, due to Einhorn)
# Line 15 | Line 16 | static char SCCSid[] = "$SunId$ LBL";
16   *      vert_ill -      Total vertical illuminance
17   *
18   *              12 April 1991   Greg Ward       EPFL
19 + *              19 April 1993   R. Compagnon    EPFL (added dgi, brs_gi, ugr)
20   */
21 +
22 + #include <string.h>
23  
24   #include "standard.h"
25   #include "view.h"
26 <
27 < extern double   erfc();
24 <
26 >
27 >
28   double  posindex();
29 <
29 >
30   double  direct(), guth_dgr(), guth_vcp(), cie_cgi(),
31 <        indirect(), total();
32 <
31 >        indirect(), total(), dgi(), brs_gi(), ugr();
32 >
33   struct named_func {
34          char    *name;
35          double  (*func)();
36          char    *descrip;
37   } all_funcs[] = {
38 +        {"dgi", dgi, "Daylight Glare Index"},
39 +        {"brs_gi", brs_gi, "BRS Glare Index"},
40 +        {"ugr", ugr, "Unified Glare Rating"},
41          {"guth_vcp", guth_vcp, "Guth Visual Comfort Probability"},
42          {"cie_cgi", cie_cgi, "CIE Glare Index (Einhorn)"},
43          {"guth_dgr", guth_dgr, "Guth Disability Glare Rating"},
# Line 40 | Line 46 | struct named_func {
46          {"vert_ind", indirect, "Indirect Vertical Illuminance"},
47          {NULL}
48   };
49 <
49 >
50   struct glare_src {
51          FVECT   dir;            /* source direction */
52          double  dom;            /* solid angle */
53          double  lum;            /* average luminance */
54          struct glare_src        *next;
55   } *all_srcs = NULL;
56 <
56 >
57   struct glare_dir {
58          double  ang;            /* angle (in radians) */
59          double  indirect;       /* indirect illuminance */
60          struct glare_dir        *next;
61   } *all_dirs = NULL;
62 <
62 >
63   #define newp(type)      (type *)malloc(sizeof(type))
64 <
64 >
65   char    *progname;
66   int     print_header = 1;
67 <
67 >
68   VIEW    midview = STDVIEW;
69 <
69 >
70   int     wrongformat = 0;
71 <
72 <
71 >
72 >
73   main(argc, argv)
74   int     argc;
75   char    *argv[];
76   {
71 –        extern char     *rindex();
77          struct named_func       *funp;
78          char    *progtail;
79          int     i;
80                                          /* get program name */
81          progname = argv[0];
82 <        progtail = rindex(progname, '/');       /* final component */
82 >        progtail = strrchr(progname, '/');      /* final component */
83          if (progtail == NULL)
84                  progtail = progname;
85          else
# Line 118 | Line 123 | userr:
123                  fprintf(stderr, "\t%12s\t%s\n", funp->name, funp->descrip);
124          exit(1);
125   }
126 <
127 <
126 >
127 >
128 > int
129   headline(s)                     /* get line from header */
130   char    *s;
131   {
132          char    fmt[32];
133 <
133 >
134          if (print_header)               /* copy to output */
135                  fputs(s, stdout);
136          if (isview(s))
# Line 133 | Line 139 | char   *s;
139                  formatval(fmt, s);
140                  wrongformat = strcmp(fmt, "ascii");
141          }
142 +        return(0);
143   }
144 <
145 <
144 >
145 >
146   init()                          /* initialize calculation */
147   {
148                                          /* read header */
# Line 151 | Line 158 | init()                         /* initialize calculation */
158                  exit(1);
159          }
160   }
161 <
162 <
161 >
162 >
163   read_input()                    /* read glare sources from stdin */
164   {
165   #define S_SEARCH        0
# Line 162 | Line 169 | read_input()                   /* read glare sources from stdin */
169          char    buf[128];
170          register struct glare_src       *gs;
171          register struct glare_dir       *gd;
172 <
172 >
173          while (fgets(buf, sizeof(buf), stdin) != NULL)
174                  switch (state) {
175                  case S_SEARCH:
# Line 212 | Line 219 | readerr:
219   #undef S_SOURCE
220   #undef S_DIREC
221   }
222 <
223 <
222 >
223 >
224   print_values(funp)              /* print out calculations */
225   double  (*funp)();
226   {
227          register struct glare_dir       *gd;
228 <
228 >
229          for (gd = all_dirs; gd != NULL; gd = gd->next)
230                  printf("%f\t%f\n", gd->ang*(180.0/PI), (*funp)(gd));
231   }
232 <
233 <
232 >
233 >
234   double
235   direct(gd)                      /* compute direct vertical illuminance */
236   struct glare_dir        *gd;
# Line 231 | Line 238 | struct glare_dir       *gd;
238          FVECT   mydir;
239          double  d, dval;
240          register struct glare_src       *gs;
241 <
241 >
242          spinvector(mydir, midview.vdir, midview.vup, gd->ang);
243          dval = 0.0;
244          for (gs = all_srcs; gs != NULL; gs = gs->next) {
# Line 241 | Line 248 | struct glare_dir       *gd;
248          }
249          return(dval);
250   }
251 <
252 <
251 >
252 >
253   double
254   indirect(gd)                    /* return indirect vertical illuminance */
255   struct glare_dir        *gd;
256   {
257          return(gd->indirect);
258   }
259 <
260 <
259 >
260 >
261   double
262   total(gd)                       /* return total vertical illuminance */
263   struct glare_dir        *gd;
264   {
265          return(direct(gd)+gd->indirect);
266   }
267 <
268 <
267 >
268 >
269   /*
270   * posindex -   compute glare position index from:
271   *
# Line 272 | Line 279 | struct glare_dir       *gd;
279   * This calculation presumes the view direction and up vectors perpendicular.
280   * We return a value less than zero for improper positions.
281   */
282 <
282 >
283   double
284   posindex(sd, vd, vu)                    /* compute position index */
285   FVECT   sd, vd, vu;
286   {
287          double  sigma, tau;
288          double  d;
289 <
289 >
290          d = DOT(sd,vd);
291          if (d <= 0.0)
292                  return(-1.0);
293          if (d >= 1.0)
294                  return(1.0);
295          sigma = acos(d) * (180./PI);
296 <        d = DOT(sd,vu)/sqrt(1.0-d*d);
296 >        d = fabs(DOT(sd,vu)/sqrt(1.0-d*d));
297          if (d >= 1.0)
298                  tau = 0.0;
299          else
# Line 295 | Line 302 | FVECT  sd, vd, vu;
302                          + sigma*(21. + tau*(.26667 + tau*-.002963))*1e-5 )
303                  ) );
304   }
305 +
306 +
307 + double
308 + dgi(gd)         /* compute Daylight Glare Index */
309 + struct glare_dir        *gd;
310 + {
311 +        register struct glare_src       *gs;
312 +        FVECT   mydir,testdir[7],vhor;
313 +        double  r,posn,omega,p[7],sum;
314 +        int     i,n;
315 +
316 +        spinvector(mydir, midview.vdir, midview.vup, gd->ang);
317 +        sum = 0.0; n = 0;
318 +        for (gs = all_srcs; gs != NULL; gs = gs->next) {
319  
320 +                /* compute 1/p^2 weighted solid angle of the source */
321 +                r = sqrt(1 - pow(1.-gs->dom/2./PI,2.));
322 +                fcross(vhor,gs->dir,midview.vup);
323 +                normalize(vhor);
324 +                VCOPY(testdir[0],gs->dir);
325 +                fvsum(testdir[1],gs->dir,vhor,r);
326 +                fvsum(testdir[2],gs->dir,vhor,0.5*r);
327 +                fvsum(testdir[5],testdir[2],midview.vup,-0.866*r);
328 +                fvsum(testdir[2],testdir[2],midview.vup,0.866*r);
329 +                fvsum(testdir[3],gs->dir,vhor,-r);
330 +                fvsum(testdir[4],gs->dir,vhor,-0.5*r);
331 +                fvsum(testdir[6],testdir[4],midview.vup,0.866*r);
332 +                fvsum(testdir[4],testdir[4],midview.vup,-0.866*r);
333 +                for (i = 0; i < 7; i++) {
334 +                        normalize(testdir[i]);
335 +                        posn = posindex(testdir[i],mydir,midview.vup);
336 +                        if (posn <= FTINY)
337 +                                p[i] = 0.0;
338 +                        else
339 +                                p[i] = 1./(posn*posn);
340 +                }
341 +                r = 1-gs->dom/2./PI;
342 +                omega = gs->dom*p[0];
343 +                omega += (r*PI*(1+1/r/r)-2*PI)*(-p[0]+(p[1]+p[2])*0.5);
344 +                omega += (2*PI-r*PI*(1+1/r/r))*(-p[0]-0.1667*(p[1]+p[3])
345 +                          +0.3334*(p[2]+p[4]+p[5]+p[6]));
346  
347 +                sum += pow(gs->lum,1.6) * pow(omega,0.8) /
348 +                       (gd->indirect/PI + 0.07*sqrt(gs->dom)*gs->lum);
349 +                n++;
350 +        }
351 +        if (n == 0)
352 +                return(0.0);
353 +        return( 10*log10(0.478*sum) );
354 + }
355 +
356 +
357   double
358 + brs_gi(gd)              /* compute BRS Glare Index */
359 + struct glare_dir        *gd;
360 + {
361 +        register struct glare_src       *gs;
362 +        FVECT   mydir;
363 +        double  p;
364 +        double  sum;
365 +
366 +        spinvector(mydir, midview.vdir, midview.vup, gd->ang);
367 +        sum = 0.0;
368 +        for (gs = all_srcs; gs != NULL; gs = gs->next) {
369 +                p = posindex(gs->dir, mydir, midview.vup);
370 +                if (p <= FTINY)
371 +                        continue;
372 +                sum += pow(gs->lum/p,1.6) * pow(gs->dom,0.8);
373 +        }
374 +        if (sum <= FTINY)
375 +                return(0.0);
376 +        sum /= gd->indirect/PI;
377 +        return(10*log10(0.478*sum));
378 + }
379 +
380 +
381 + double
382   guth_dgr(gd)            /* compute Guth discomfort glare rating */
383   struct glare_dir        *gd;
384   {
# Line 308 | Line 389 | struct glare_dir       *gd;
389          double  sum;
390          double  wtot, brsum;
391          int     n;
392 <
392 >
393          spinvector(mydir, midview.vdir, midview.vup, gd->ang);
394          sum = wtot = brsum = 0.0; n = 0;
395          for (gs = all_srcs; gs != NULL; gs = gs->next) {
# Line 322 | Line 403 | struct glare_dir       *gd;
403          }
404          if (n == 0)
405                  return(0.0);
325 –
406          return( pow(.5*sum/pow((brsum+(5.-wtot)*gd->indirect/PI)/5.,.44),
407                          pow((double)n, -.0914) ) );
408   #undef q
409   }
410 <
411 <
332 < extern double   erf(), erfc();
333 <
410 >
411 >
412   #ifndef M_SQRT2
413   #define M_SQRT2 1.41421356237309504880
414   #endif
415 <
415 >
416   #define norm_integral(z)        (1.-.5*erfc((z)/M_SQRT2))
417 <
418 <
417 >
418 >
419   double
420   guth_vcp(gd)            /* compute Guth visual comfort probability */
421   struct glare_dir        *gd;
422   {
423 +        extern double   erfc();
424          double  dgr;
425 <
425 >
426          dgr = guth_dgr(gd);
427          if (dgr <= FTINY)
428                  return(100.0);
429          return(100.*norm_integral(6.374-1.3227*log(dgr)));
430   }
431 <
432 <
431 >
432 >
433   double
434   cie_cgi(gd)             /* compute CIE Glare Index */
435   struct glare_dir        *gd;
# Line 360 | Line 439 | struct glare_dir       *gd;
439          double  dillum;
440          double  p;
441          double  sum;
442 <
442 >
443          spinvector(mydir, midview.vdir, midview.vup, gd->ang);
444          sum = 0.0;
445          for (gs = all_srcs; gs != NULL; gs = gs->next) {
# Line 373 | Line 452 | struct glare_dir       *gd;
452                  return(0.0);
453          dillum = direct(gd);
454          return(8.*log10(2.*sum*(1.+dillum/500.)/(dillum+gd->indirect)));
455 + }
456 +
457 +
458 + double
459 + ugr(gd)         /* compute Unified Glare Rating */
460 + struct glare_dir        *gd;
461 + {
462 +        register struct glare_src       *gs;
463 +        FVECT   mydir;
464 +        double  p;
465 +        double  sum;
466 +
467 +        spinvector(mydir, midview.vdir, midview.vup, gd->ang);
468 +        sum = 0.0;
469 +        for (gs = all_srcs; gs != NULL; gs = gs->next) {
470 +                p = posindex(gs->dir, mydir, midview.vup);
471 +                if (p <= FTINY)
472 +                        continue;
473 +                sum += gs->lum*gs->lum * gs->dom / (p*p);
474 +        }
475 +        if (sum <= FTINY)
476 +                return(0.0);
477 +        return(8.*log10(0.25*sum*PI/gd->indirect));
478   }

Diff Legend

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