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

Comparing ray/src/util/evalglare.c (file contents):
Revision 2.4 by greg, Thu Jul 28 16:46:10 2016 UTC vs.
Revision 2.9 by greg, Fri Nov 23 19:32:17 2018 UTC

# Line 1 | Line 1
1   #ifndef lint
2   static const char RCSid[] = "$Id$";
3   #endif
4 < /* EVALGLARE V1.29
4 > /* EVALGLARE V2.07
5   * Evalglare Software License, Version 2.0
6   *
7   * Copyright (c) 1995 - 2016 Fraunhofer ISE, EPFL.
# Line 310 | Line 310 | changed masking threshold to 0.05 cd/m2
310     */
311   /* evalglare.c, v1.30 2016/07/28 change zonal output: If just one zone is specified, only one zone shows up in the output (bug removal).
312     */
313 + /* evalglare.c, v1.31 2016/08/02  bug removal: default output did not calculate the amout of glare sources before and therefore no_glaresources was set to zero causing dgi,ugr being set to zero as well. Now renumbering of the glare sources and calculation of the amount of glare sources is done for all output versions.
314 +   */
315 + /* evalglare.c, v2.00 2016/11/15  add of a second fast calculation mode for annual calculations, activted by -2. Output: dgp,ugr
316 +   */
317 + /* evalglare.c, v2.01 2016/11/16  change of -2 option (now -2 dir_illum). External provision of the direct illuminance necessary, since the sun interpolation of daysim is causing problems in calculation of the background luminance.
318 +   */
319 + /* evalglare.c, v2.02 2017/02/28  change of warning message, when invalid exposure setting is found. Reason: tab removal is not in all cases the right measure - it depends which tool caused the invalid exposure entry   */
320 +
321 + /* evalglare.c, v2.03 2017/08/04  ad of -O option - disk replacement by providing luminance, not documented
322 +  */
323 +
324 + /* evalglare.c, v2.04 2017/08/04  adding -q option: use of Ev/pi as background luminance. not documented. no combination with -n option!!!
325 +  */
326 +
327 + /* evalglare.c, v2.05 2018/08/28  change of the -q option for the choice of the background luminance calculation mode: 0: CIE method (Ev-Edir)/pi, 1: mathematical correct average background luminance, 2: Ev/pi.
328 + change of default options:
329 + - cosinus weighted calculation of the background luminance (according to CIE) is now default.
330 + - absolute threshold for the glare source detection is now default (2000cd/m2), based on study of C. Pierson
331 +  */
332 +
333 + /* evalglare.c, v2.06 2018/08/29  
334 + change of default value of multiplier b to 5.0, if task options (-t or -T ) are activated AND -b NOT used. To be downward compatible when using the task method.
335 +  */
336 +
337 + /* evalglare.c, v2.07 2018/11/06  
338 + bugfix: correction of error in the equations of PGSV_con and PGSV_sat
339 + all three PGSV equations are calculated now
340 + illuminance from the masking area (E_v_mask) is also printed
341 + bugfix: in VCPs error fuction equation, value of 6.347 replaced by 6.374
342 +  */
343 +  
344   #define EVALGLARE
345   #define PROGNAME "evalglare"
346 < #define VERSION "1.30 release 29.07.2016 by EPFL, J.Wienold"
346 > #define VERSION "2.07 release 17.11.2018 by EPFL, J.Wienold"
347   #define RELEASENAME PROGNAME " " VERSION
348  
349  
# Line 1164 | Line 1195 | float get_vcp(float dgr )
1195   {
1196          float vcp;
1197  
1198 <        vcp = 50 * erf((6.347 - 1.3227 * log(dgr)) / 1.414213562373) + 50;
1198 >        vcp = 50 * erf((6.374 - 1.3227 * log(dgr)) / 1.414213562373) + 50;
1199          if (dgr > 750) {
1200                  vcp = 0;
1201          }
# Line 1199 | Line 1230 | float get_ugr(pict * p, double lum_backg, int igs, int
1230                  }
1231          }
1232          ugr = 8 * log10(0.25 / lum_backg * sum_glare);
1233 <
1233 >        if (sum_glare==0) {
1234 >        ugr=0.0;
1235 >        }
1236 >        if (lum_backg<=0) {
1237 >        ugr=-99.0;
1238 >        }
1239 >        
1240          return ugr;
1241  
1242   }
# Line 1329 | Line 1366 | float get_cgi(pict * p, double E_v, double E_v_dir, in
1366          return cgi;
1367   }      
1368  
1369 < /* subroutine for the calculation of the PGSV; is only applied, when masking is done, since it refers only to the window. Important: masking area must be the window! */
1370 < float get_pgsv(double E_v, double E_mask,double omega_mask,double lum_mask_av)
1369 >
1370 >
1371 > /* subroutine for the calculation of the PGSV_con; is only applied, when masking is done, since it refers only to the window. Important: masking area must be the window! */
1372 > float get_pgsv_con(double E_v, double E_mask,double omega_mask,double lum_mask_av, double Lavg)
1373   {
1374 <        float pgsv;
1374 >        float pgsv_con;
1375          double Lb;
1376  
1377 <        Lb = (E_v-E_mask)/1.414213562373;
1377 > /*        Lb = (E_v-E_mask)/3.14159265359;  */
1378 > /*        Lb = (2*E_v-lum_mask_av*omega_mask)/(2*3.14159265359-omega_mask); */
1379 >          Lb = (2*3.14159265359*Lavg-lum_mask_av*omega_mask)/(2*3.14159265359-omega_mask);
1380  
1340 –        pgsv =3.2*log10(lum_mask_av)-0.64*log10(omega_mask)+(0.79*log10(omega_mask)-0.61)*log10(Lb)-8.2 ;
1381  
1382 +        pgsv_con =3.2*log10(lum_mask_av)-0.64*log10(omega_mask)+(0.79*log10(omega_mask)-0.61)*log10(Lb)-8.2 ;
1383  
1384 <        return pgsv;
1384 >
1385 >        return pgsv_con;
1386   }
1387  
1388   /* subroutine for the calculation of the PGSV_saturation; is only applied, when masking is done, since it refers only to the window. Important: masking area must be the window! */
# Line 1348 | Line 1390 | float get_pgsv_sat(double E_v)
1390   {
1391          float pgsv_sat;
1392  
1393 <        pgsv_sat =3.3-(0.57+3.3)/pow((1+E_v/1.414213562373/1250),1.7);
1393 >        pgsv_sat =3.3-(0.57+3.3)/(1+pow(E_v/3.14159265359/1250,1.7));
1394  
1395  
1396          return pgsv_sat;
1397   }
1398  
1399 + /* subroutine for the calculation of the PGSV; is only applied, when masking is done, since it refers only to the window. Important: masking area must be the window! */
1400  
1401 + float get_pgsv(double E_v, double E_mask,double omega_mask,double lum_mask_av,double Ltask, double Lavg)
1402 + {
1403 +        float pgsv;
1404 +        double Lb;
1405  
1406 + /*        Lb = (E_v-E_mask)/3.14159265359;  */
1407 + /*        Lb = (2*E_v-lum_mask_av*omega_mask)/(2*3.14159265359-omega_mask); */
1408 +          Lb = (2*3.14159265359*Lavg-lum_mask_av*omega_mask)/(2*3.14159265359-omega_mask);
1409 +        
1410 +        if (Lb==0.0 ) {
1411 +               fprintf(stderr,  " warning: Background luminance is 0 or masking area = full image! pgsv cannot be calculated (set to -99)!!\n");
1412 +                pgsv=-99;
1413 +                        }else{
1414 +                if ( (lum_mask_av/Lb) > (E_v/(3.14159265359*Ltask))) {
1415 +                        pgsv=get_pgsv_con(E_v,E_mask,omega_mask,lum_mask_av, Lavg);
1416 +                }else{
1417 +                        pgsv=get_pgsv_sat(E_v)  ;
1418 +                        }}
1419 +        return pgsv;
1420  
1421 + }
1422 +
1423 +
1424 +
1425   #ifdef  EVALGLARE
1426  
1427  
# Line 1371 | Line 1436 | int main(int argc, char **argv)
1436          pict *pm = pict_create();
1437          int     skip_second_scan,calcfast,age_corr,cut_view,cut_view_type,calc_vill, output, detail_out2, x1,y1, fill, yfillmax, yfillmin,
1438                  ext_vill, set_lum_max, set_lum_max2, img_corr,x_disk,y_disk,task_color, i_splitstart,zones,act_gsn,splitgs,
1439 <                i_split, posindex_2, task_lum, checkfile, rval, i, i_max, x, y,x2,y2,x_zone,y_zone, i_z1, i_z2,
1439 >                i_split, posindex_2, task_lum, checkfile, rval, i, i_max, x, y,x2,y2,x_zone,y_zone, i_z1, i_z2, thres_activate,
1440                  igs, actual_igs, lastpixelwas_gs, icol, xt, yt, change,checkpixels, before_igs, sgs, splithigh,uniform_gs,x_max, y_max,y_mid,
1441                  detail_out, posindex_picture, non_cos_lb, rx, ry, rmx,rmy,apply_disability,band_calc,band_color,masking,i_mask,no_glaresources,force;
1442 <        double  lum_total_max,age_corr_factor,age,dgp_ext,dgp,low_light_corr,omega_cos_contr, setvalue, lum_ideal, E_v_contr, sigma,om,delta_E,
1443 <                E_vl_ext, lum_max, new_lum_max, r_center, ugp, ugr_exp, dgi_mod,lum_a, pgsv,E_v_mask,pgsv_sat,angle_disk,dist,n_corner_px,zero_corner_px,
1444 <                search_pix, a1, a2, a3, a4, a5, c3, c1, c2, r_split, max_angle,r_actual,lum_actual,
1442 >        double  LUM_replace,lum_total_max,age_corr_factor,age,dgp_ext,dgp,low_light_corr,omega_cos_contr, setvalue, lum_ideal, E_v_contr, sigma,om,delta_E,
1443 >                E_vl_ext, lum_max, new_lum_max, r_center, ugp, ugr_exp, dgi_mod,lum_a, E_v_mask,angle_disk,dist,n_corner_px,zero_corner_px,
1444 >                search_pix, a1, a2, a3, a4, a5, c3, c1, c2, r_split, max_angle,r_actual,lum_actual,dir_ill,
1445                  omegat, sang, E_v, E_v2, E_v_dir, avlum, act_lum, ang, angle_z1, angle_z2,per_95_band,per_75_band,pos,
1446                  l_max, lum_backg, lum_backg_cos, omega_sources, lum_sources,per_75_mask,per_95_mask,per_75_z1,per_95_z1,per_75_z2,per_95_z2,
1447                  lum, lum_source,teta,Lveil_cie,Lveil_cie_sum,disability_thresh,u_r,u_g,u_b,band_angle,band_avlum,
# Line 1385 | Line 1450 | int main(int argc, char **argv)
1450                  lum_z1[1],lum_z1_av,lum_z1_std[1],lum_z1_median[1],omega_z1,bbox_z1[2],
1451                  lum_z2[1],lum_z2_av,lum_z2_std[1],lum_z2_median[1],omega_z2,bbox_z2[2],
1452                  lum_pos[1],lum_nopos_median[1],lum_pos_median[1],lum_pos2_median[1],lum_pos_mean,lum_pos2_mean;
1453 <        float lum_task, lum_thres, dgi,  vcp, cgi, ugr, limit, dgr,
1453 >        float lum_task, lum_thres, dgi,  vcp, cgi, ugr, limit, dgr,pgsv ,pgsv_sat,pgsv_con,
1454                  abs_max, Lveil;
1455          char maskfile[500],file_out[500], file_out2[500], version[500];
1456          char *cline;
# Line 1425 | Line 1490 | int main(int argc, char **argv)
1490          muc_rvar_clear(s_posweight2);
1491  
1492          /*set required user view parameters to invalid values*/
1493 +        dir_ill=0.0;
1494          delta_E=0.0;
1495          no_glaresources=0;
1496          n_corner_px=0;
# Line 1446 | Line 1512 | int main(int argc, char **argv)
1512          lum_band_av = 0.0;
1513          omega_band = 0.0;
1514          pgsv = 0.0 ;
1515 +        pgsv_con = 0.0 ;
1516          pgsv_sat = 0.0 ;
1517          E_v_mask = 0.0;
1518          lum_z1_av = 0.0;
# Line 1487 | Line 1554 | int main(int argc, char **argv)
1554          omega_cos_contr = 0.0;
1555          lum_ideal = 0.0;
1556          max_angle = 0.2;
1557 <        lum_thres = 5.0;
1557 >        lum_thres = 2000.0;
1558 >        lum_task = 0.0;
1559          task_lum = 0;
1560          sgs = 0;
1561          splithigh = 1;
# Line 1505 | Line 1573 | int main(int argc, char **argv)
1573          c1 = 5.87e-05;
1574          c2 = 0.092;
1575          c3 = 0.159;
1576 <        non_cos_lb = 1;
1576 >        non_cos_lb = 0;
1577          posindex_2 = 0;
1578          task_color = 0;
1579          limit = 50000.0;
# Line 1527 | Line 1595 | int main(int argc, char **argv)
1595          omega_mask=0.0;
1596          i_mask=0;
1597          actual_igs=0;
1598 +        LUM_replace=0;
1599 +        thres_activate=0;
1600   /* command line for output picture*/
1601  
1602          cline = (char *) malloc(CLINEMAX+1);
# Line 1577 | Line 1647 | int main(int argc, char **argv)
1647                          break;
1648                  case 'b':
1649                          lum_thres = atof(argv[++i]);
1650 +                        thres_activate = 1;
1651                          break;
1652                  case 'c':
1653                          checkfile = 1;
# Line 1665 | Line 1736 | int main(int argc, char **argv)
1736                          strcpy(file_out2, argv[++i]);
1737   /*                      printf("max lum set to %f\n",new_lum_max);*/
1738                          break;
1739 +                case 'O':
1740 +                        img_corr = 1;
1741 +                        set_lum_max2 = 3;
1742 +                        x_disk = atoi(argv[++i]);
1743 +                        y_disk = atoi(argv[++i]);
1744 +                        angle_disk = atof(argv[++i]);
1745 +                        LUM_replace = atof(argv[++i]);
1746 +                        strcpy(file_out2, argv[++i]);
1747 + /*                      printf("max lum set to %f\n",new_lum_max);*/
1748 +                        break;
1749  
1750 <                case 'n':
1750 >
1751 > /* deactivated          case 'n':
1752                          non_cos_lb = 0;
1753                          break;
1754 + */
1755 +                case 'q':
1756 +                        non_cos_lb = atoi(argv[++i]);
1757 +                        break;
1758  
1759                  case 't':
1760                          task_lum = 1;
# Line 1737 | Line 1823 | int main(int argc, char **argv)
1823                  case '1':
1824                          output = 1;
1825                          break;
1826 +                case '2':
1827 +                        output = 2;
1828 +                        dir_ill = atof(argv[++i]);
1829 +                        break;
1830  
1831                  case 'v':
1832                          if (argv[i][2] == '\0') {
# Line 1766 | Line 1856 | int main(int argc, char **argv)
1856                  }
1857          }
1858  
1859 + /* set multiplier for task method to 5, if not specified */
1860 +
1861 + if ( task_lum == 1 && thres_activate == 0){
1862 +                lum_thres = 5.0;
1863 + }
1864   /*fast calculation, if gendgp_profile is used: No Vertical illuminance calculation, only dgp is calculated*/
1865  
1866 < if (output == 1 && ext_vill == 1) {
1866 > if (output == 1 && ext_vill == 1 ) {
1867                         calcfast=1;
1868                         }
1869 +                      
1870 + if (output == 2 && ext_vill == 1 ) {
1871 +                       calcfast=2;
1872 +                       }
1873 +                      
1874   /*masking and zoning cannot be applied at the same time*/
1875  
1876   if (masking ==1 && zones >0) {
# Line 2059 | Line 2159 | if (cut_view==2) {
2159          lum_pos_mean= lum_pos_mean/sang;
2160          lum_pos2_mean= lum_pos2_mean/sang;
2161  
2162 <        if (set_lum_max2 >= 1 && E_v_contr > 0 && (E_vl_ext - E_v) > 0) {
2162 >        if ((set_lum_max2 >= 1 && E_v_contr > 0 && (E_vl_ext - E_v) > 0 ) || set_lum_max2==3) {
2163  
2164 +                if (set_lum_max2<3){
2165                  lum_ideal = (E_vl_ext - E_v + E_v_contr) / omega_cos_contr;
2166                  if (set_lum_max2 == 2 && lum_ideal >= 2e9) {
2167                  printf("warning! luminance of replacement pixels would be larger than 2e9 cd/m2. Value set to 2e9cd/m2!\n") ;
# Line 2071 | Line 2172 | if (cut_view==2) {
2172                  }
2173                  printf("change luminance values!! lum_ideal,setvalue,E_vl_ext,E_v,E_v_contr %f  %f %f %f %f\n",
2174                           lum_ideal, setvalue, E_vl_ext, E_v, E_v_contr);
2175 +                          }else{setvalue=LUM_replace;
2176 +                         }
2177  
2178              
2179                  for (x = 0; x < pict_get_xsize(p); x++)
# Line 2090 | Line 2193 | if (cut_view==2) {
2193                                                          pict_get_color(p, x, y)[BLU] =
2194                                                                  setvalue / 179.0;
2195  
2196 <                                                }else{ if(set_lum_max2 ==2 ) {
2196 >                                                }else{ if(set_lum_max2 >1 ) {
2197                                                          r_actual =acos(DOT(pict_get_cached_dir(p, x_disk, y_disk), pict_get_cached_dir(p, x, y))) * 2;
2198                                                          if (x_disk == x && y_disk==y ) r_actual=0.0;
2199  
# Line 2103 | Line 2206 | if (cut_view==2) {
2206                                                          pict_get_color(p, x, y)[BLU] =
2207                                                                  setvalue / 179.0;
2208                                                        
2209 <                                                       }
2107 <                                                                
2108 <                                                
2209 >                                                       }                                                
2210                                                  }
2211                                                  }
2212                                          }
# Line 2221 | Line 2322 | if (cut_view==1) {
2322   /*      pict_write(p,"firstscan.pic");   */
2323  
2324  
2325 < if (calcfast == 1 || search_pix <= 1.0) {
2325 >
2326 >
2327 > if (calcfast ==1 || search_pix <= 1.0 || calcfast == 2 ) {
2328     skip_second_scan=1;
2329     }
2330 +  
2331  
2332   /* second glare source scan: combine glare sources facing each other */
2333          change = 1;
# Line 2398 | Line 2502 | if (calcfast == 1 || search_pix <= 1.0) {
2502  
2503   /* calculation of direct vertical illuminance for CGI and for disability glare, coloring glare sources*/
2504  
2505 <        if (calcfast == 0) {
2505 >        if (calcfast == 0 || calcfast == 2) {
2506          for (x = 0; x < pict_get_xsize(p); x++)
2507                  for (y = 0; y < pict_get_ysize(p); y++) {
2508                          if (pict_get_hangle(p, x, y, p->view.vdir, p->view.vup, &ang)) {
# Line 2486 | Line 2590 | if (calcfast == 1 || search_pix <= 1.0) {
2590                          lum_backg = lum_backg_cos;
2591          }
2592  
2593 +        if (non_cos_lb == 2) {
2594 +                        lum_backg = E_v / 3.1415927;
2595 +        }
2596 +
2597 +
2598   /* file writing NOT here
2599          if (checkfile == 1) {
2600                  pict_write(p, file_out);
# Line 2495 | Line 2604 | if (calcfast == 1 || search_pix <= 1.0) {
2604   /* print detailed output */
2605          
2606          
2498 –        if (detail_out == 1) {
2607  
2608   /* masking */
2609  
# Line 2538 | Line 2646 | if (calcfast == 1 || search_pix <= 1.0) {
2646                  muc_rvar_get_median(s_mask,lum_mask_median);
2647                  muc_rvar_get_bounding_box(s_mask,bbox);
2648   /* PSGV only why masking of window is applied! */
2649 <                 pgsv = get_pgsv(E_v, E_v_mask, omega_mask, lum_mask_av);
2649 >
2650 >        
2651 >        if (task_lum == 0 || lum_task == 0.0 ) {
2652 >                        fprintf(stderr, " warning: Task area not set or task luminance=0 ! pgsv cannot be calculated (set to -99)!!\n");
2653 >                        pgsv = -99;
2654 >                } else {
2655 >                        pgsv = get_pgsv(E_v, E_v_mask, omega_mask, lum_mask_av,lum_task,avlum);
2656 >                        }
2657 >
2658 >                 pgsv_con = get_pgsv_con(E_v, E_v_mask, omega_mask, lum_mask_av,avlum);
2659                   pgsv_sat =get_pgsv_sat(E_v);
2543 –                printf ("masking:no_pixels,omega,av_lum,median_lum,std_lum,perc_75,perc_95,lum_min,lum_max,pgsv,pgsv_sat: %i %f %f %f %f %f %f %f %f %f %f\n",i_mask,omega_mask,lum_mask_av,lum_mask_median[0],sqrt(lum_mask_std[0]),per_75_mask,per_95_mask,bbox[0],bbox[1],pgsv,pgsv_sat );
2660  
2661 +        if (detail_out == 1) {
2662 +
2663 +                printf ("masking:no_pixels,omega,av_lum,median_lum,std_lum,perc_75,perc_95,lum_min,lum_max,pgsv_con,pgsv_sat,pgsv,Ev_mask: %i %f %f %f %f %f %f %f %f %f %f %f %f\n",i_mask,omega_mask,lum_mask_av,lum_mask_median[0],sqrt(lum_mask_std[0]),per_75_mask,per_95_mask,bbox[0],bbox[1],pgsv_con,pgsv_sat,pgsv,E_v_mask );
2664 +
2665 +        }      
2666                  
2546 –                
2667          }
2668  
2669   /* zones */
# Line 2633 | Line 2753 | if (calcfast == 1 || search_pix <= 1.0) {
2753                  per_95_z1=lum_z1_median[0];
2754                  muc_rvar_get_median(s_z1,lum_z1_median);
2755                  muc_rvar_get_bounding_box(s_z1,bbox_z1);
2756 <
2756 >        if (detail_out == 1) {
2757 >
2758                  printf ("zoning:z1_omega,z1_av_lum,z1_median_lum,z1_std_lum,z1_perc_75,z1_perc_95,z1_lum_min,z1_lum_max: %f %f %f %f %f %f %f %f\n",omega_z1,lum_z1_av,lum_z1_median[0],sqrt(lum_z1_std[0]),per_75_z1,per_95_z1,bbox_z1[0],bbox_z1[1] );
2759  
2760                 if (zones == 2 ) {
2761  
2762                  printf ("zoning:z2_omega,z2_av_lum,z2_median_lum,z2_std_lum,z2_perc_75,z2_perc_95,z2_lum_min,z2_lum_max:  %f %f %f %f %f %f %f %f\n",omega_z2,lum_z2_av,lum_z2_median[0],sqrt(lum_z2_std[0]),per_75_z2,per_95_z2,bbox_z2[0],bbox_z2[1] );
2763 < }            
2763 > } }            
2764                  
2765          }
2766  
# Line 2653 | Line 2774 | if (calcfast == 1 || search_pix <= 1.0) {
2774                  no_glaresources=i;
2775  
2776   /* glare sources */
2777 +        if (detail_out == 1) {
2778 +
2779                  printf
2780                          ("%i No pixels x-pos y-pos L_s Omega_s Posindx L_b L_t E_vert Edir Max_Lum Sigma xdir ydir zdir Eglare_cie Lveil_cie teta glare_zone\n",
2781                           i);
# Line 2817 | Line 2940 | if (calcfast == 0) {
2940                                  if (E_vl_ext < 1000) {
2941                                  low_light_corr=1.0*exp(0.024*E_vl_ext-4)/(1+exp(0.024*E_vl_ext-4)); } else {low_light_corr=1.0 ;}
2942                                  dgp =low_light_corr*dgp;
2943 <                                dgp =age_corr_factor*dgp;
2944 <                printf("%f\n", dgp);
2943 >
2944 >                     if (calcfast == 2) {
2945 >                    
2946 >                         lum_backg_cos=(E_vl_ext-dir_ill)/3.1415927;
2947 >                         ugr = get_ugr(p, lum_backg_cos, igs, posindex_2);
2948 >                         printf("%f %f \n", dgp,ugr);
2949 >                     }else{      
2950 >                         printf("%f\n", dgp);
2951 >                }
2952          }
2953  
2954  

Diff Legend

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