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

Comparing ray/src/util/evalglare.c (file contents):
Revision 2.6 by greg, Thu May 18 02:25:27 2017 UTC vs.
Revision 3.06 by greg, Mon Sep 15 20:21:24 2025 UTC

# Line 1 | Line 1
1   #ifndef lint
2 < static const char RCSid[] = "$Id$";
2 > static const char RCSid[] = "$Id$";
3   #endif
4 < /* EVALGLARE V2.00
5 < * Evalglare Software License, Version 2.0
4 > /* EVALGLARE V3.05
5 > * Evalglare Software License, Version 3.0x
6   *
7 < * Copyright (c) 1995 - 2016 Fraunhofer ISE, EPFL.
7 > * Copyright (c) 1995 - 2022 Fraunhofer ISE, EPFL.
8   * All rights reserved.
9   *
10   *
# Line 317 | Line 317 | changed masking threshold to 0.05 cd/m2
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 <  
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/17  
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 > /* evalglare.c, v2.08 2018/11/27
345 > bugfix: checkroutine for same image size for the masking corrected
346 >  */
347 >
348 > /* evalglare.c, v2.09 2019/01/18
349 > calculate "illuminance-contribution of zones"
350 > -switch to turn off low-light correction: 4
351 >  */
352 >
353 > /* evalglare.c, v2.10 2019/07/30, 2020/06/25
354 > -add multiimage-mode for annual glare evaluation -Q
355 > -extra output for multiimage mode
356 > 2020/06/22 - added Ev contribution from each glare source as output to the multiimage mode
357 > -switch for correction modes (-C), value of 0 switches off all corrections
358 >  */
359 >
360 >
361 > /* evalglare.c, v2.11  2020/07/06
362 > -add extra output for DGP-contribution in -d mode
363 > -fix (remove) 1 pixel offset of the mulltiimage-mode
364 >  */
365 >
366 > /* evalglare.c, v2.12  2021/07/19
367 > -bugfix: center of glare source calculation now based on l*omega weighting, before it was just on omega.
368 >  */
369 >
370 > /* evalglare.c, v2.13  2021/07/21
371 > - center of glare source calculation now based on l*l*omega weighting.
372 > - bugfix: Position-index below line of sight is now corrected and based on the modified equation from Iwata and Osterhaus 2010
373 >  */
374 >
375 > /* evalglare.c, v2.14  2021/10/20
376 > - center of glare source calculation now based on l*omega weighting.
377 > - adding UGP2 equation to the output
378 > - new default setting: low-light correction off !
379 >  */
380 >
381 > /* evalglare.c, v3.0  2021/12/24
382 > - adding patch-mode (-z patchangle combineflag ) as glare source grouping method
383 >  */
384 >
385 > /* evalglare.c, v3.01  2022/01/05
386 > - change of multi-image-mode. Adding extra scans with scaled image
387 > syntax now: -Q n_images n_extrascans_per_image n_images*( imagename x y  n_extrascans*(scale) )
388 >  */
389 > /* evalglare.c, v3.02  2022/02/11
390 > - correction of position index below line of sight after Takuro Kikuchi's comments in radiance-discourse from Dec 2021.  
391 >  */
392 >
393 > /* evalglare.c, v3.03  2022/03/18, 2024/06/25
394 > - add no pixels output for zonal analysis.
395 >
396 > /* evalglare.c, v3.04  2023/09/01
397 > - test version to include peripherical dependency, changes not (yet) used for following versions
398 >
399 > /* evalglare.c, v3.05  2024/06/25
400 > - fix bug of still having lowlight correction even if turned off
401 > - changing abs_max,Lveil from float to double
402 >  */
403 >
404 >
405   #define EVALGLARE
406   #define PROGNAME "evalglare"
407 < #define VERSION "2.02 release 28.02.2017 by EPFL, J.Wienold"
407 > #define VERSION "3.05 release 25.06.2024 by J.Wienold, EPFL"
408   #define RELEASENAME PROGNAME " " VERSION
409  
410  
# Line 331 | Line 415 | changed masking threshold to 0.05 cd/m2
415   #include "platform.h"
416   #include "muc_randvar.h"
417  
334 char *progname;
335
418   /* subroutine to add a pixel to a glare source */
419   void add_pixel_to_gs(pict * p, int x, int y, int gsn)
420   {
421          double old_av_posx, old_av_posy, old_av_lum, old_omega, act_omega,
422 <                new_omega, act_lum;
422 >                new_omega, act_lum,temp_av_posx,temp_av_posy,av_lum;
423  
424  
425          pict_get_npix(p, gsn) = pict_get_npix(p, gsn) + 1;
# Line 349 | Line 431 | void add_pixel_to_gs(pict * p, int x, int y, int gsn)
431          act_omega = pict_get_omega(p, x, y);
432          act_lum = luminance(pict_get_color(p, x, y));
433          new_omega = old_omega + act_omega;
434 <        pict_get_av_posx(p, gsn) =
435 <                (old_av_posx * old_omega + x * act_omega) / new_omega;
436 <        pict_get_av_posy(p, gsn) =
437 <                (old_av_posy * old_omega + y * act_omega) / new_omega;
434 >        pict_get_av_lum(p, gsn) =
435 >                (old_av_lum * old_omega + act_lum * act_omega) / new_omega;
436 >
437 >        av_lum=pict_get_av_lum(p, gsn);
438 >        temp_av_posx =
439 >                (old_av_posx *  old_omega* old_av_lum + x * act_omega*act_lum) / (old_av_lum*old_omega + act_lum* act_omega );
440 >        pict_get_av_posx(p, gsn) = temp_av_posx;
441 >        temp_av_posy =
442 >                (old_av_posy *  old_omega* old_av_lum + y * act_omega*act_lum) / (old_av_lum*old_omega + act_lum* act_omega );
443 >        
444 >        pict_get_av_posy(p, gsn) = temp_av_posy;
445          if (isnan((pict_get_av_posx(p, gsn))))
446                  fprintf(stderr,"error in add_pixel_to_gs %d %d %f %f %f %f\n", x, y, old_av_posy, old_omega, act_omega, new_omega);
447  
359        pict_get_av_lum(p, gsn) =
360                (old_av_lum * old_omega + act_lum * act_omega) / new_omega;
448          pict_get_av_omega(p, gsn) = new_omega;
449          pict_get_gsn(p, x, y) = gsn;
450          if (act_lum < pict_get_lum_min(p, gsn)) {
# Line 731 | Line 818 | void add_secondary_gs(pict * p, int x, int y, float r,
818   void split_pixel_from_gs(pict * p, int x, int y, int new_gsn, int uniform_gs, double u_r, double u_g , double u_b)
819   {
820          int old_gsn, icol;
821 <        double old_av_posx, old_av_posy, old_av_lum, old_omega, act_omega,
822 <                new_omega, act_lum;
821 >        double old_av_posx, old_av_posy, old_av_lum, old_omega, act_omega, delta_E,
822 >                new_omega, act_lum,temp_av_posx,temp_av_posy,av_lum;
823  
824  
825   /* change existing gs */
# Line 743 | Line 830 | void split_pixel_from_gs(pict * p, int x, int y, int n
830          old_av_posx = pict_get_av_posx(p, old_gsn);
831          old_av_posy = pict_get_av_posy(p, old_gsn);
832          old_omega = pict_get_av_omega(p, old_gsn);
833 +        
834 +        
835 +        
836  
837          new_omega = old_omega - act_omega;
838          pict_get_av_omega(p, old_gsn) = new_omega;
749        pict_get_av_posx(p, old_gsn) =
750                (old_av_posx * old_omega - x * act_omega) / new_omega;
751        pict_get_av_posy(p, old_gsn) =
752                (old_av_posy * old_omega - y * act_omega) / new_omega;
839  
754
840          act_lum = luminance(pict_get_color(p, x, y));
841          old_av_lum = pict_get_av_lum(p, old_gsn);
842          pict_get_av_lum(p, old_gsn) =
843                  (old_av_lum * old_omega - act_lum * act_omega) / new_omega;
844 +
845 +        av_lum = pict_get_av_lum(p, old_gsn);
846 +        temp_av_posx =
847 +                (old_av_posx *old_av_lum* old_omega  - x *act_lum* act_omega ) / (old_av_lum*old_omega  - act_lum* act_omega);
848 +
849 +        pict_get_av_posx(p, old_gsn) = temp_av_posx;
850 +        temp_av_posy =
851 +                (old_av_posy *old_av_lum* old_omega - y *act_lum* act_omega ) / (old_av_lum*old_omega  - act_lum* act_omega);
852 +        pict_get_av_posy(p, old_gsn) = temp_av_posy;
853 +
854          /* add pixel to new  gs */
855  
856          add_pixel_to_gs(p, x, y, new_gsn);
857  
858 +
859 +
860 +
861   /* color pixel of new gs */
862  
863   /*      icol = setglcolor(p, x, y, new_gsn, uniform_gs, u_r, u_g , u_b); */
# Line 771 | Line 869 | void split_pixel_from_gs(pict * p, int x, int y, int n
869   float get_posindex(pict * p, float x, float y, int postype)
870   {
871          float posindex;
872 <        double teta, phi, sigma, tau, deg, d, s, r, fact;
872 >        double teta, beta, phi, sigma, tau, deg, d, s, r, fact;
873  
874  
875          pict_get_vangle(p, x, y, p->view.vdir, p->view.vup, &phi);
# Line 799 | Line 897 | float get_posindex(pict * p, float x, float y, int pos
897   /* KIM  model */        
898          posindex = exp ((sigma-(-0.000009*tau*tau*tau+0.0014*tau*tau+0.0866*tau+21.633))/(-0.000009*tau*tau*tau+0.0013*tau*tau+0.0853*tau+8.772));
899          }else{
900 +
901   /* Guth model, equation from IES lighting handbook */
902          posindex =
903                  exp((35.2 - 0.31889 * tau -
# Line 807 | Line 906 | float get_posindex(pict * p, float x, float y, int pos
906                                                                                                                   0.002963 * tau *
907                                                                                                                   tau) / 100000 *
908                          sigma * sigma);
810 /* below line of sight, using Iwata model */
811        if (phi < 0) {
812                d = 1 / tan(phi);
813                s = tan(teta) / tan(phi);
814                r = sqrt(1 / d * 1 / d + s * s / d / d);
815                if (r > 0.6)
816                        fact = 1.2;
817                if (r > 3) {
818                        fact = 1.2;
819                        r = 3;
820                }
909  
910 <                posindex = 1 + fact * r;
910 > /* below line of sight, using Iwata model, CIE2010, converted coordinate system according to Takuro Kikuchi */
911 >
912 >        if (phi < 0) {
913 >          beta = atan(tan(sigma/deg)* sqrt(1 + 0.3225 * pow(cos(tau/deg),2))) * deg;
914 >          posindex = exp(6.49 / 1000 * beta + 21.0 / 100000 * beta * beta);
915 >
916          }
917 +
918                  if (posindex > 16)
919                  posindex = 16;
920   }
# Line 1114 | Line 1208 | get_dgp(pict * p, double E_v, int igs, double a1, doub
1208                                                                                             pict_get_av_posy(p, i),
1209                                                                                             posindex_2),
1210                                                                    a4) * pow(pict_get_av_omega(p, i), a2);
1211 < /*                      printf("i,sum_glare %i %f\n",i,sum_glare);*/                                      
1211 >                        /*printf("i,sum_glare %i %f\n",i,sum_glare);    */                                
1212                          }
1213                  }
1214                  dgp =
# Line 1172 | Line 1266 | float get_vcp(float dgr )
1266   {
1267          float vcp;
1268  
1269 <        vcp = 50 * erf((6.347 - 1.3227 * log(dgr)) / 1.414213562373) + 50;
1269 >        vcp = 50 * erf((6.374 - 1.3227 * log(dgr)) / 1.414213562373) + 50;
1270          if (dgr > 750) {
1271                  vcp = 0;
1272          }
# Line 1268 | Line 1362 | float get_ugp(pict * p, double lum_backg, int igs, int
1362          return ugp;
1363  
1364   }
1365 + /* subroutine for the calculation of the updated unified glare probability (Hirning 2016)*/
1366 + float get_ugp2(pict * p, double lum_backg, int igs, int posindex_2)
1367 + {
1368 +        float ugp;
1369 +        double sum_glare,ugp2;
1370 +        int i, i_glare;
1371  
1372  
1373 +        sum_glare = 0;
1374 +        i_glare = 0;
1375 +        for (i = 0; i <= igs; i++) {
1376 +
1377 +                if (pict_get_npix(p, i) > 0) {
1378 +                        i_glare = i_glare + 1;
1379 +                        sum_glare +=
1380 +                                pow(pict_get_av_lum(p, i) /
1381 +                                        get_posindex(p, pict_get_av_posx(p, i),
1382 +                                                                 pict_get_av_posy(p, i), posindex_2),
1383 +                                        2) * pict_get_av_omega(p, i);
1384 +
1385 +                }
1386 +        }
1387 +        ugp2 = 1/pow(1.0+2.0/7.0*pow(sum_glare/lum_backg,-0.2),10.0);
1388 +
1389 +        return ugp2;
1390 +
1391 + }
1392 +
1393 +
1394   /* subroutine for the calculation of the disability glare according to Poynter */
1395   float get_disability(pict * p, double lum_backg, int igs)
1396   {
# Line 1343 | Line 1464 | float get_cgi(pict * p, double E_v, double E_v_dir, in
1464          return cgi;
1465   }      
1466  
1467 < /* 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! */
1468 < float get_pgsv(double E_v, double E_mask,double omega_mask,double lum_mask_av)
1467 >
1468 >
1469 > /* 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! */
1470 > float get_pgsv_con(double E_v, double E_mask,double omega_mask,double lum_mask_av, double Lavg)
1471   {
1472 <        float pgsv;
1472 >        float pgsv_con;
1473          double Lb;
1474  
1475 <        Lb = (E_v-E_mask)/1.414213562373;
1475 > /*        Lb = (E_v-E_mask)/3.14159265359;  */
1476 > /*        Lb = (2*E_v-lum_mask_av*omega_mask)/(2*3.14159265359-omega_mask); */
1477 >          Lb = (2*3.14159265359*Lavg-lum_mask_av*omega_mask)/(2*3.14159265359-omega_mask);
1478  
1354        pgsv =3.2*log10(lum_mask_av)-0.64*log10(omega_mask)+(0.79*log10(omega_mask)-0.61)*log10(Lb)-8.2 ;
1479  
1480 +        pgsv_con =3.2*log10(lum_mask_av)-0.64*log10(omega_mask)+(0.79*log10(omega_mask)-0.61)*log10(Lb)-8.2 ;
1481  
1482 <        return pgsv;
1482 >
1483 >        return pgsv_con;
1484   }
1485  
1486   /* 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 1362 | Line 1488 | float get_pgsv_sat(double E_v)
1488   {
1489          float pgsv_sat;
1490  
1491 <        pgsv_sat =3.3-(0.57+3.3)/pow((1+E_v/1.414213562373/1250),1.7);
1491 >        pgsv_sat =3.3-(0.57+3.3)/(1+pow(E_v/3.14159265359/1250,1.7));
1492  
1493  
1494          return pgsv_sat;
1495   }
1496  
1497 + /* 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! */
1498  
1499 + float get_pgsv(double E_v, double E_mask,double omega_mask,double lum_mask_av,double Ltask, double Lavg)
1500 + {
1501 +        float pgsv;
1502 +        double Lb;
1503  
1504 + /*        Lb = (E_v-E_mask)/3.14159265359;  */
1505 + /*        Lb = (2*E_v-lum_mask_av*omega_mask)/(2*3.14159265359-omega_mask); */
1506 +          Lb = (2*3.14159265359*Lavg-lum_mask_av*omega_mask)/(2*3.14159265359-omega_mask);
1507 +        
1508 +        if (Lb==0.0 ) {
1509 +               fprintf(stderr,  " warning: Background luminance is 0 or masking area = full image! pgsv cannot be calculated (set to -99)!!\n");
1510 +                pgsv=-99;
1511 +                        }else{
1512 +                if ( (lum_mask_av/Lb) > (E_v/(3.14159265359*Ltask))) {
1513 +                        pgsv=get_pgsv_con(E_v,E_mask,omega_mask,lum_mask_av, Lavg);
1514 +                }else{
1515 +                        pgsv=get_pgsv_sat(E_v)  ;
1516 +                        }}
1517 +        return pgsv;
1518  
1519 + }
1520 +
1521 +
1522 +
1523   #ifdef  EVALGLARE
1524  
1525  
# Line 1380 | Line 1529 | float get_pgsv_sat(double E_v)
1529  
1530   int main(int argc, char **argv)
1531   {
1532 <        #define CLINEMAX 4095 /* memory allocated for command line string */
1532 >        #define CLINEMAX 999999999 /* memory allocated for command line string */
1533          pict *p = pict_create();
1534          pict *pm = pict_create();
1535 <        int     skip_second_scan,calcfast,age_corr,cut_view,cut_view_type,calc_vill, output, detail_out2, x1,y1, fill, yfillmax, yfillmin,
1536 <                ext_vill, set_lum_max, set_lum_max2, img_corr,x_disk,y_disk,task_color, i_splitstart,zones,act_gsn,splitgs,
1537 <                i_split, posindex_2, task_lum, checkfile, rval, i, i_max, x, y,x2,y2,x_zone,y_zone, i_z1, i_z2,
1535 >        pict *pcoeff = pict_create();
1536 >        int     add,i1,i2,y_lines,ix_offset,x_offset,y_offset, patchmode,ix, iy, patch_pixdistance_x, patch_pixdistance_y,nx_patch,ny_patch, lowlight,skip_second_scan,
1537 >                calcfast,age_corr,cut_view,cut_view_type,calc_vill, output, detail_out2, x1,y1, fill, yfillmax, yfillmin,
1538 >                ext_vill, set_lum_max, set_lum_max2, img_corr,x_disk,y_disk,task_color, i_splitstart,zones,act_gsn,splitgs,j,jj,multi_image_mode,
1539 >                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,num_images,xmap,ymap,
1540                  igs, actual_igs, lastpixelwas_gs, icol, xt, yt, change,checkpixels, before_igs, sgs, splithigh,uniform_gs,x_max, y_max,y_mid,
1541 <                detail_out, posindex_picture, non_cos_lb, rx, ry, rmx,rmy,apply_disability,band_calc,band_color,masking,i_mask,no_glaresources,force;
1542 <        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,
1543 <                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,
1541 >                detail_out, posindex_picture, non_cos_lb, rx, ry,lastpixelwas_peak, rmx,rmy,apply_disability,band_calc,band_color,masking,i_mask,no_glaresources,force,num_scans;
1542 >        double  abs_max,Lveil,xxx,angle_v,angle_h,patch_angle, r_contrast,r_dgp,r_glare,sum_glare, LUM_replace,lum_total_max,age_corr_factor,age,dgp_ext,dgp,low_light_corr,omega_cos_contr, setvalue,
1543 >                lum_ideal, E_v_contr, sigma,om,delta_E,
1544 >                E_vl_ext, lum_max, new_lum_max, r_center, ugp,ugp2, ugr_exp, dgi_mod,lum_a, E_v_mask,angle_disk,dist,n_corner_px,zero_corner_px,
1545                  search_pix, a1, a2, a3, a4, a5, c3, c1, c2, r_split, max_angle,r_actual,lum_actual,dir_ill,
1546 <                omegat, sang, E_v, E_v2, E_v_dir, avlum, act_lum, ang, angle_z1, angle_z2,per_95_band,per_75_band,pos,
1546 >                omegat, sang,Ez1,Ez2, E_v, E_v2, E_v_dir, avlum, act_lum, ang, angle_z1, angle_z2,per_95_band,per_75_band,pos,
1547                  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,
1548                  lum, lum_source,teta,Lveil_cie,Lveil_cie_sum,disability_thresh,u_r,u_g,u_b,band_angle,band_avlum,
1549                  lum_mask[1],lum_mask_av,lum_mask_std[1],lum_mask_median[1],omega_mask,bbox[2],
# Line 1399 | Line 1551 | int main(int argc, char **argv)
1551                  lum_z1[1],lum_z1_av,lum_z1_std[1],lum_z1_median[1],omega_z1,bbox_z1[2],
1552                  lum_z2[1],lum_z2_av,lum_z2_std[1],lum_z2_median[1],omega_z2,bbox_z2[2],
1553                  lum_pos[1],lum_nopos_median[1],lum_pos_median[1],lum_pos2_median[1],lum_pos_mean,lum_pos2_mean;
1554 <        float lum_task, lum_thres, dgi,  vcp, cgi, ugr, limit, dgr,
1555 <                abs_max, Lveil;
1404 <        char maskfile[500],file_out[500], file_out2[500], version[500];
1554 >        float lum_task, lum_thres, dgi,  vcp, cgi, ugr, limit, dgr,pgsv ,pgsv_sat,pgsv_con;
1555 >        char maskfile[500],file_out[500], file_out2[500], version[500],correction_type[500];
1556          char *cline;
1557 +        char** temp_image_name;
1558 +        int *x_temp_img, *y_temp_img;
1559 +        double *scale_image_scans;
1560          VIEW userview = STDVIEW;
1561          int gotuserview = 0;
1562 +        
1563          struct muc_rvar* s_mask;
1564          s_mask = muc_rvar_create();
1565          muc_rvar_set_dim(s_mask, 1);
# Line 1439 | Line 1594 | int main(int argc, char **argv)
1594          muc_rvar_clear(s_posweight2);
1595  
1596          /*set required user view parameters to invalid values*/
1597 <        dir_ill=0.0;
1597 >        patchmode=0;
1598 >        patch_angle=25.0;
1599 >        patch_pixdistance_x=0;
1600 >        patch_pixdistance_y=0;
1601 >        lowlight=0;
1602 >        multi_image_mode=0;
1603 >        lastpixelwas_peak=0;    
1604 >        num_images=0;        
1605 >        dir_ill=0.0;
1606          delta_E=0.0;
1607          no_glaresources=0;
1608          n_corner_px=0;
# Line 1461 | Line 1624 | int main(int argc, char **argv)
1624          lum_band_av = 0.0;
1625          omega_band = 0.0;
1626          pgsv = 0.0 ;
1627 +        pgsv_con = 0.0 ;
1628          pgsv_sat = 0.0 ;
1629          E_v_mask = 0.0;
1630 +        Ez1 = 0.0;
1631 +        Ez2 = 0.0;
1632          lum_z1_av = 0.0;
1633          omega_z1 = 0.0;
1634          lum_z2_av = 0.0;
# Line 1502 | Line 1668 | int main(int argc, char **argv)
1668          omega_cos_contr = 0.0;
1669          lum_ideal = 0.0;
1670          max_angle = 0.2;
1671 <        lum_thres = 5.0;
1671 >        lum_thres = 2000.0;
1672 >        lum_task = 0.0;
1673          task_lum = 0;
1674          sgs = 0;
1675          splithigh = 1;
# Line 1520 | Line 1687 | int main(int argc, char **argv)
1687          c1 = 5.87e-05;
1688          c2 = 0.092;
1689          c3 = 0.159;
1690 <        non_cos_lb = 1;
1690 >        non_cos_lb = 0;
1691          posindex_2 = 0;
1692          task_color = 0;
1693          limit = 50000.0;
# Line 1528 | Line 1695 | int main(int argc, char **argv)
1695          set_lum_max2 = 0;
1696          img_corr=0;
1697          abs_max = 0;
1698 <        progname = argv[0];
1698 >        fixargv0(argv[0]);
1699          E_v_contr = 0.0;
1700 <        strcpy(version, "1.19 release 09.12.2015 by J.Wienold");
1700 >        strcpy(version, "3.05 release 25.06.2024 by J.Wienold");
1701          low_light_corr=1.0;
1702          output = 0;
1703          calc_vill = 0;
# Line 1542 | Line 1709 | int main(int argc, char **argv)
1709          omega_mask=0.0;
1710          i_mask=0;
1711          actual_igs=0;
1712 +        LUM_replace=0;
1713 +        thres_activate=0;
1714   /* command line for output picture*/
1715  
1716 <        cline = (char *) malloc(CLINEMAX+1);
1716 >        cline = (char *) malloc(CLINEMAX+100);
1717          cline[0] = '\0';
1718          for (i = 0; i < argc; i++) {
1719   /*       fprintf(stderr, "%d %d \n",i,strlen(argv[i]));*/
# Line 1592 | Line 1761 | int main(int argc, char **argv)
1761                          break;
1762                  case 'b':
1763                          lum_thres = atof(argv[++i]);
1764 +                        lum_source =lum_thres;
1765 +                        thres_activate = 1;
1766                          break;
1767                  case 'c':
1768                          checkfile = 1;
# Line 1606 | Line 1777 | int main(int argc, char **argv)
1777                  case 'r':
1778                          max_angle = atof(argv[++i]);
1779                          break;
1780 +                case 'z':
1781 +                        patch_angle = atof(argv[++i]);
1782 +                        patchmode= atoi(argv[++i]);
1783 +                        if ( patchmode == 3) { output=4;}
1784 +                        
1785 +                        /* patchmode = 0 : deactivated; patchmode =1: patch without combining, normal output; patchmode =3: patch without combining, coefficient output; patchmode =2 patch with combining (not yet implemented...) */
1786 +                        break;
1787 +
1788                  case 's':
1789                          sgs = 1;
1790                          break;
# Line 1680 | Line 1859 | int main(int argc, char **argv)
1859                          strcpy(file_out2, argv[++i]);
1860   /*                      printf("max lum set to %f\n",new_lum_max);*/
1861                          break;
1862 +                case 'O':
1863 +                        img_corr = 1;
1864 +                        set_lum_max2 = 3;
1865 +                        x_disk = atoi(argv[++i]);
1866 +                        y_disk = atoi(argv[++i]);
1867 +                        angle_disk = atof(argv[++i]);
1868 +                        LUM_replace = atof(argv[++i]);
1869 +                        strcpy(file_out2, argv[++i]);
1870 + /*                      printf("max lum set to %f\n",new_lum_max);*/
1871 +                        break;
1872  
1873 <                case 'n':
1873 >
1874 > /* deactivated          case 'n':
1875                          non_cos_lb = 0;
1876                          break;
1877 + */
1878 +                case 'q':
1879 +                        non_cos_lb = atoi(argv[++i]);
1880 +                        break;
1881  
1882                  case 't':
1883                          task_lum = 1;
# Line 1748 | Line 1942 | int main(int argc, char **argv)
1942                          /*case 'v':
1943                          printf("evalglare  %s \n",version);
1944                          exit(1); */
1945 +                case 'C':
1946 +                        strcpy(correction_type,argv[++i]);
1947 +                        
1948 +                        if (!strcmp(correction_type,"l-")  ){
1949 +                /*      printf("low light off!\n"); */
1950 +                                                                                        lowlight = 0; }
1951 +                        if (!strcmp(correction_type,"l+")  ){
1952 +                /*      printf("low light on!\n"); */
1953 +                                                                                        lowlight = 1; }
1954 +                        if (!strcmp(correction_type,"0")  ){
1955 +                /*      printf("all corrections off!\n"); */
1956 +                                                                                        lowlight = 0; }
1957 +                        
1958 +                        break;
1959  
1960 +                        /*case 'v':
1961 +                        printf("evalglare  %s \n",version);
1962 +                        exit(1); */
1963 +
1964                  case '1':
1965                          output = 1;
1966                          break;
# Line 1756 | Line 1968 | int main(int argc, char **argv)
1968                          output = 2;
1969                          dir_ill = atof(argv[++i]);
1970                          break;
1971 <
1971 >                case '3':
1972 >                        output = 3;
1973 >                        break;
1974 >                case '4':
1975 >                        lowlight = 0;
1976 >                        break;
1977 >                case 'Q':
1978 >                        multi_image_mode=1;
1979 >                        output= 3;                      
1980 >                        calcfast=1;
1981 >                        num_images =atoi(argv[++i]);
1982 >                        num_scans =atoi(argv[++i]);
1983 >                        temp_image_name = malloc(sizeof(char*)*num_images);
1984 >                        
1985 >                        x_temp_img=(int *) malloc(sizeof(int) * num_images);
1986 >                        y_temp_img=(int *) malloc(sizeof(int) * num_images);
1987 >                        scale_image_scans=(double *) malloc(sizeof(double) * (num_scans+1)*num_images);
1988 >        /* iterate through all images and allocate 256 characters to each: */
1989 >                        for (j = 0; j < num_images; j++) {
1990 >                                scale_image_scans[j*(num_scans+1)]=1.0;
1991 >                                temp_image_name[j] = malloc(256*sizeof(char));
1992 >                                strcpy(temp_image_name[j], argv[++i]);
1993 >                                x_temp_img[j] = atoi(argv[++i]);
1994 >                                y_temp_img[j] = atoi(argv[++i]);
1995 >                                for (jj=1;jj<=num_scans;jj++){
1996 >                                      scale_image_scans[j*(num_scans+1)+jj]=atof(argv[++i]);}
1997 >                        }
1998 >                
1999 >                
2000 >                        break;
2001                  case 'v':
2002                          if (argv[i][2] == '\0') {
2003                               printf("%s \n",RELEASENAME);                              
# Line 1785 | Line 2026 | int main(int argc, char **argv)
2026                  }
2027          }
2028  
2029 + /* set multiplier for task method to 5, if not specified */
2030 +                
2031 +
2032 +
2033 + if ( task_lum == 1 && thres_activate == 0){
2034 +                lum_thres = 5.0;
2035 + }
2036   /*fast calculation, if gendgp_profile is used: No Vertical illuminance calculation, only dgp is calculated*/
2037  
2038   if (output == 1 && ext_vill == 1 ) {
# Line 1849 | Line 2097 | if (masking ==1 && zones >0) {
2097                  exit(1);
2098          }
2099  
2100 +        if ( patchmode > 0 && p->view.type != VT_ANG) {
2101 +        
2102 +                fprintf(stderr, "error: Patchmode only possible with view type vta  !  Stopping... \n");
2103 +                exit(1);
2104 +        
2105 +        }
2106 +        
2107 +        
2108          if (p->view.type == VT_PER) {
2109                  fprintf(stderr, "warning: image has perspective view type specified in header ! \n");
2110          }
2111  
2112          if (masking == 1) {
2113  
2114 <                if (!pict_get_xsize(p)==pict_get_xsize(pm) || !pict_get_ysize(p)==pict_get_ysize(pm)) {
2114 >                if (pict_get_xsize(p)!=pict_get_xsize(pm) || pict_get_ysize(p)!=pict_get_ysize(pm)) {
2115                  fprintf(stderr, "error: masking image has other resolution than main image ! \n");
2116                  fprintf(stderr, "size must be identical \n");
2117                  printf("resolution main image : %dx%d\n",pict_get_xsize(p),pict_get_ysize(p));
# Line 1890 | Line 2146 | if (masking ==1 && zones >0) {
2146  
2147                  }
2148          }
2149 <
2149 >        
2150 >
2151   /* Check task position  */
2152  
2153          if (task_lum == 1) {
# Line 1914 | Line 2171 | if (masking ==1 && zones >0) {
2171          pict_get_z1_gsn(p,igs) = 0;
2172          pict_get_z2_gsn(p,igs) = 0;
2173  
2174 + if (multi_image_mode<1) {
2175  
2176 +
2177   /* cut out GUTH field of view and exit without glare evaluation */
2178   if (cut_view==2) {
2179          if (cut_view_type==1) {
# Line 2083 | Line 2342 | if (cut_view==2) {
2342          lum_pos_mean= lum_pos_mean/sang;
2343          lum_pos2_mean= lum_pos2_mean/sang;
2344  
2345 <        if (set_lum_max2 >= 1 && E_v_contr > 0 && (E_vl_ext - E_v) > 0) {
2345 >        if ((set_lum_max2 >= 1 && E_v_contr > 0 && (E_vl_ext - E_v) > 0 ) || set_lum_max2==3) {
2346  
2347 +                if (set_lum_max2<3){
2348                  lum_ideal = (E_vl_ext - E_v + E_v_contr) / omega_cos_contr;
2349                  if (set_lum_max2 == 2 && lum_ideal >= 2e9) {
2350                  printf("warning! luminance of replacement pixels would be larger than 2e9 cd/m2. Value set to 2e9cd/m2!\n") ;
# Line 2095 | Line 2355 | if (cut_view==2) {
2355                  }
2356                  printf("change luminance values!! lum_ideal,setvalue,E_vl_ext,E_v,E_v_contr %f  %f %f %f %f\n",
2357                           lum_ideal, setvalue, E_vl_ext, E_v, E_v_contr);
2358 +                          }else{setvalue=LUM_replace;
2359 +                         }
2360  
2361              
2362                  for (x = 0; x < pict_get_xsize(p); x++)
# Line 2114 | Line 2376 | if (cut_view==2) {
2376                                                          pict_get_color(p, x, y)[BLU] =
2377                                                                  setvalue / 179.0;
2378  
2379 <                                                }else{ if(set_lum_max2 ==2 ) {
2379 >                                                }else{ if(set_lum_max2 >1 ) {
2380                                                          r_actual =acos(DOT(pict_get_cached_dir(p, x_disk, y_disk), pict_get_cached_dir(p, x, y))) * 2;
2381                                                          if (x_disk == x && y_disk==y ) r_actual=0.0;
2382  
# Line 2127 | Line 2389 | if (cut_view==2) {
2389                                                          pict_get_color(p, x, y)[BLU] =
2390                                                                  setvalue / 179.0;
2391                                                        
2392 <                                                       }
2131 <                                                                
2132 <                                                
2392 >                                                       }                                                
2393                                                  }
2394                                                  }
2395                                          }
# Line 2201 | Line 2461 | if (cut_view==1) {
2461   /*     printf("Ev, avlum, Omega_tot%f %f %f \n",E_v, avlum, sang); */
2462  
2463   /* first glare source scan: find primary glare sources */
2464 +
2465 + if (patchmode==0) {
2466          for (x = 0; x < pict_get_xsize(p); x++) {
2467 < /*                lastpixelwas_gs=0; */
2467 >                lastpixelwas_gs=0;
2468 > /*              lastpixelwas_peak=0; */
2469                  for (y = 0; y < pict_get_ysize(p); y++) {
2470                          if (pict_get_hangle(p, x, y, p->view.vdir, p->view.vup, &ang)) {
2471                                  if (pict_is_validpixel(p, x, y)) {
# Line 2212 | Line 2475 | if (cut_view==1) {
2475                                                                               lum_total_max=act_lum;
2476                                                                                 }
2477   /* speed improvement first scan: when last pixel was glare source, then use this glare source number instead of searching.
2478 <   has been deactivated with v1.25
2478 >   has been deactivated with v1.25, reactivated with v2.10 */
2479                        
2480 <                                                if (lastpixelwas_gs==0 || search_pix <= 1.0 ) { */
2480 >                                                if (lastpixelwas_gs==0 || search_pix <= 1.0 ) {
2481                                                  actual_igs = find_near_pgs(p, x, y, max_angle, 0, igs, 1);
2482 < /*  }*/
2482 >  }
2483                                                  if (actual_igs == 0) {
2484                                                          igs = igs + 1;
2485                                                          pict_new_gli(p);
# Line 2232 | Line 2495 | if (cut_view==1) {
2495                                                  pict_get_gsn(p, x, y) = actual_igs;
2496                                                  pict_get_pgs(p, x, y) = 1;
2497                                                  add_pixel_to_gs(p, x, y, actual_igs);
2498 < /*                                                lastpixelwas_gs=actual_igs; */
2498 >                                                lastpixelwas_gs=actual_igs;
2499  
2500                                          } else {
2501                                                  pict_get_pgs(p, x, y) = 0;
2502 < /*                                               lastpixelwas_gs=0;*/
2502 >                                                lastpixelwas_gs=0;
2503                                          }
2504                                  }
2505                          }
2506                  }
2507               }
2508 < /*      pict_write(p,"firstscan.pic");   */
2508 >         } else {    
2509 > /* patchmode on!
2510 > calculation only for angular projection!
2511  
2512 + */
2513  
2514  
2515 + angle_v=p->view.vert;
2516 + angle_h=p->view.horiz;
2517  
2518 < if (calcfast ==1 || search_pix <= 1.0 || calcfast == 2 ) {
2518 >
2519 > /*printf ("angle_v,angle_h: %f %f \n",angle_v,angle_h) ;
2520 >
2521 >  setting up patches as glare sources */
2522 >
2523 > patch_pixdistance_x=floor(pict_get_xsize(p)/angle_h*patch_angle);
2524 > patch_pixdistance_y=floor(pict_get_ysize(p)/angle_v*patch_angle);
2525 >
2526 > nx_patch=floor(angle_v/patch_angle)+1;
2527 > ny_patch=floor(angle_h/patch_angle)+1;
2528 >
2529 > y_offset=floor (patch_pixdistance_y/2);
2530 > x_offset=floor (patch_pixdistance_x/2);
2531 > /* printf ("nx_patch,ny_patch,x_offset,y_offset,patch_pixdistance_x,patch_pixdistance_y   %i %i %i %i %i %i\n",nx_patch,ny_patch,x_offset,y_offset,patch_pixdistance_x,patch_pixdistance_y) ; */
2532 >
2533 > ix_offset=0;
2534 > for (iy=1; iy<=ny_patch;iy++) {
2535 >
2536 >
2537 >
2538 > for (ix=1; ix<=nx_patch;ix++) {
2539 >                                                        igs = igs + 1;
2540 >                                                        pict_new_gli(p);
2541 >
2542 >                                                        pict_get_lum_min(p, igs) = HUGE_VAL;
2543 >                                                        pict_get_Eglare(p,igs) = 0.0;
2544 >                                                        pict_get_Dglare(p,igs) = 0.0;
2545 >                                                        pict_get_z1_gsn(p,igs) = 0;
2546 >                                                        pict_get_z2_gsn(p,igs) = 0;
2547 >                                                        pict_get_dx_max(p,igs) = (x_offset+ix_offset*x_offset+(ix-1)*patch_pixdistance_x)*1.0;
2548 >                                                        pict_get_dy_max(p,igs) = (y_offset+(iy-1)*patch_pixdistance_y)*1.0;
2549 >
2550 > /* printf ("igs, x-patch, y-patch : %i %f %f \n",igs,pict_get_dx_max(p,igs),pict_get_dy_max(p,igs) ) ; */
2551 >
2552 > }
2553 > ix_offset=ix_offset+1;
2554 > if (ix_offset==2) {
2555 >                        ix_offset =0 ;
2556 >                        }
2557 >
2558 > }
2559 >                for (y = 0; y < pict_get_ysize(p); y++) {
2560 >                        for (x = 0; x < pict_get_xsize(p); x++) {
2561 >
2562 >                        if (pict_get_hangle(p, x, y, p->view.vdir, p->view.vup, &ang)) {
2563 >                                if (pict_is_validpixel(p, x, y)) {
2564 >                                        act_lum = luminance(pict_get_color(p, x, y));
2565 >                                        if (act_lum > lum_source) {
2566 >                                                if (act_lum >lum_total_max) {
2567 >                                                                             lum_total_max=act_lum;
2568 >                                                                               }
2569 >                                        
2570 >                                                y_lines = floor((y)/patch_pixdistance_y);
2571 >                                                xxx = (x+0.0)/(patch_pixdistance_x+0.0)-0.5*(y_lines % 2);
2572 >                                                if (xxx<0 ) { xxx=0.0 ;}
2573 >                                                i1 = y_lines*(nx_patch)+floor(xxx)+1;
2574 >                                                i2=0;
2575 >                                                add=0;
2576 >                                                if (y_lines % 2 == 1 ) {add=1;}
2577 >                                                
2578 >                                                if (y >pict_get_dy_max(p,i1)) {
2579 >                                                        
2580 >                                                        if (x > pict_get_dx_max(p,i1)) {
2581 >                                                                i2=i1+nx_patch+add;
2582 >                                                                }else {
2583 >                                                                i2=i1+nx_patch-1+add;
2584 >                                                                }
2585 >                                                        }else {
2586 >                                                
2587 >                                                        if (x > pict_get_dx_max(p,i1)) {
2588 >                                                                i2=i1-nx_patch+add;
2589 >                                                                }else {
2590 >                                                                i2=i1-nx_patch-1+add;
2591 >                                                                }
2592 >                                                        }
2593 >                                                
2594 >                                                
2595 >        
2596 >                                                
2597 >                                                if (i2 > igs || i2 < 1) {actual_igs=i1;}else{
2598 >                                                if ( ((x-pict_get_dx_max(p,i1))*(x-pict_get_dx_max(p,i1))+(y-pict_get_dy_max(p,i1))*(y-pict_get_dy_max(p,i1))) < ((x-pict_get_dx_max(p,i2))*(x-pict_get_dx_max(p,i2))+(y-pict_get_dy_max(p,i2))*(y-pict_get_dy_max(p,i2))) ) {
2599 >                                                 actual_igs=i1; }else {actual_igs=i2;}}
2600 >                                                
2601 >        
2602 >                                                pict_get_gsn(p, x, y) = actual_igs;
2603 >                                                pict_get_pgs(p, x, y) = 1;
2604 >                                                add_pixel_to_gs(p, x, y, actual_igs);
2605 >                        /*                      setglcolor(p, x, y, actual_igs, uniform_gs, u_r, u_g , u_b); */
2606 >                                }
2607 >                                }
2608 >                                }
2609 >
2610 >
2611 >
2612 >
2613 >
2614 > }               }
2615 >
2616 >
2617 >
2618 >            
2619 >             }
2620 >  /*                    pict_write(p,"firstscan.pic");   */
2621 >
2622 >
2623 >
2624 >
2625 > if (calcfast ==1 || search_pix <= 1.0 || calcfast == 2 || patchmode > 0) {
2626     skip_second_scan=1;
2627     }
2628    
# Line 2331 | Line 2706 | if (calcfast ==1 || search_pix <= 1.0 || calcfast == 2
2706  
2707   /* extract extremes from glare sources to extra glare source */
2708          if (splithigh == 1 && lum_total_max>limit) {
2709 + /*             fprintf(stderr,  " split of glare source!\n"); */
2710  
2711                  r_split = max_angle / 2.0;
2712                  for (i = 0; i <= i_max; i++) {
# Line 2490 | Line 2866 | if (calcfast ==1 || search_pix <= 1.0 || calcfast == 2
2866          i=0;
2867          for (x = 0; x <= igs; x++) {
2868          if (pict_get_npix(p, x) > 0) {
2869 +                sum_glare += pow(pict_get_av_lum(p, x),2.0)* pict_get_av_omega(p, x)/pow(get_posindex(p, pict_get_av_posx(p, x), pict_get_av_posy(p, x),posindex_2), 2.0);
2870                  lum_sources += pict_get_av_lum(p, x) * pict_get_av_omega(p, x);
2871                  omega_sources += pict_get_av_omega(p, x);
2872                  i=i+1;
2873          }}
2874 <      
2874 >        sum_glare= c2*log10(1 + sum_glare / pow(E_v, a3));
2875          if (sang == omega_sources) {
2876                 lum_backg =avlum;
2877          } else {
# Line 2513 | Line 2890 | if (calcfast ==1 || search_pix <= 1.0 || calcfast == 2
2890                          lum_backg = lum_backg_cos;
2891          }
2892  
2893 +        if (non_cos_lb == 2) {
2894 +                        lum_backg = E_v / 3.1415927;
2895 +        }
2896 +
2897 +
2898   /* file writing NOT here
2899          if (checkfile == 1) {
2900                  pict_write(p, file_out);
# Line 2564 | Line 2946 | if (calcfast ==1 || search_pix <= 1.0 || calcfast == 2
2946                  muc_rvar_get_median(s_mask,lum_mask_median);
2947                  muc_rvar_get_bounding_box(s_mask,bbox);
2948   /* PSGV only why masking of window is applied! */
2949 <                 pgsv = get_pgsv(E_v, E_v_mask, omega_mask, lum_mask_av);
2949 >
2950 >        
2951 >        if (task_lum == 0 || lum_task == 0.0 ) {
2952 >                        fprintf(stderr, " warning: Task area not set or task luminance=0 ! pgsv cannot be calculated (set to -99)!!\n");
2953 >                        pgsv = -99;
2954 >                } else {
2955 >                        pgsv = get_pgsv(E_v, E_v_mask, omega_mask, lum_mask_av,lum_task,avlum);
2956 >                        }
2957 >
2958 >                 pgsv_con = get_pgsv_con(E_v, E_v_mask, omega_mask, lum_mask_av,avlum);
2959                   pgsv_sat =get_pgsv_sat(E_v);
2960  
2961          if (detail_out == 1) {
2962  
2963 <                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 );
2963 >                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 );
2964  
2965          }      
2966                  
# Line 2597 | Line 2988 | if (calcfast ==1 || search_pix <= 1.0 || calcfast == 2
2988                                                    muc_rvar_add_sample(s_z1, lum_z1);
2989                                                    omega_z1 += pict_get_omega(p, x, y);
2990                                                    lum_z1_av += pict_get_omega(p, x, y)* lum_actual;
2991 +                                                  Ez1 += DOT(p->view.vdir, pict_get_cached_dir(p, x, y))* pict_get_omega(p, x, y)* lum_actual;
2992                                                    setglcolor(p,x,y,1,1 , 0.66, 0.01 ,0.33);
2993   /*check if separation gsn already exist */
2994  
# Line 2613 | Line 3005 | if (calcfast ==1 || search_pix <= 1.0 || calcfast == 2
3005                                                      splitgs=(int)(pict_get_z1_gsn(p,act_gsn));
3006                                          /*          printf("splitgs%i \n",splitgs);       */              
3007                                                      split_pixel_from_gs(p, x, y, splitgs, uniform_gs, u_r, u_g , u_b);
3008 +                                        /* move direct illuminance contribution into  zone -value           */
3009 +                                                delta_E=DOT(p->view.vdir, pict_get_cached_dir(p, x, y))* pict_get_omega(p, x, y)* luminance(pict_get_color(p, x, y));
3010 +                                                pict_get_Eglare(p,act_gsn ) = pict_get_Eglare(p,act_gsn ) - delta_E;
3011 +                                                pict_get_Eglare(p,igs ) = pict_get_Eglare(p,igs ) + delta_E;
3012 +            
3013 +                                                    
3014                                                  }                                
3015                                                    }
3016   /*zone2 */
# Line 2624 | Line 3022 | if (calcfast ==1 || search_pix <= 1.0 || calcfast == 2
3022                                                    muc_rvar_add_sample(s_z2, lum_z2);
3023                                                    omega_z2 +=   pict_get_omega(p, x, y);
3024                                                    lum_z2_av += pict_get_omega(p, x, y)* lum_actual;
3025 +                                                  Ez2 += DOT(p->view.vdir, pict_get_cached_dir(p, x, y))* pict_get_omega(p, x, y)* lum_actual;
3026                                                    setglcolor(p,x,y,1,1 , 0.65, 0.33 ,0.02);
3027   /*                                                printf("zone 2 x y act_gsn pict_get_z1_gsn(p,act_gsn) pict_get_z2_gsn(p,act_gsn): %i %i %f 1:%f 2:%f \n",x,y,act_gsn,pict_get_z1_gsn(p,act_gsn),pict_get_z2_gsn(p,act_gsn));
3028   */                                                if (act_gsn > 0){
# Line 2638 | Line 3037 | if (calcfast ==1 || search_pix <= 1.0 || calcfast == 2
3037                                                  splitgs=(int)(pict_get_z2_gsn(p,act_gsn));
3038   /*                                              printf("splitgs %i \n",splitgs);*/                              
3039                                                      split_pixel_from_gs(p, x, y, splitgs, uniform_gs, u_r, u_g , u_b);
3040 +                                        /* move direct illuminance contribution into  zone -value           */
3041 +                                                delta_E=DOT(p->view.vdir, pict_get_cached_dir(p, x, y))* pict_get_omega(p, x, y)* luminance(pict_get_color(p, x, y));
3042 +                                                pict_get_Eglare(p,act_gsn ) = pict_get_Eglare(p,act_gsn ) - delta_E;
3043 +                                                pict_get_Eglare(p,igs ) = pict_get_Eglare(p,igs ) + delta_E;
3044 +
3045                                             }
3046                                  }
3047  
# Line 2664 | Line 3068 | if (calcfast ==1 || search_pix <= 1.0 || calcfast == 2
3068                  muc_rvar_get_bounding_box(s_z1,bbox_z1);
3069          if (detail_out == 1) {
3070  
3071 <                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] );
3071 >                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,Ez1,i_z1: %f %f %f %f %f %f %f %f %f %i\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],Ez1,i_z1 );
3072  
3073                 if (zones == 2 ) {
3074  
3075 <                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] );
3075 >                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,Ez2,i_z2:  %f %f %f %f %f %f %f %f %f %i\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],Ez2,i_z1 );
3076   } }            
3077                  
3078          }
3079  
3080 + int new_gs_number[igs+1];
3081 +
3082                  i = 0;
3083                  for (x = 0; x <= igs; x++) {
3084   /* resorting glare source numbers */
3085                          if (pict_get_npix(p, x) > 0) {
3086                                  i = i + 1;
3087 +                                pict_get_Dglare(p,x) = i;
3088 +                                new_gs_number[x] = i;
3089 + /*                      printf("%i %i %f %i \n", i,x,pict_get_Dglare(p,x),new_gs_number[x] ); */
3090                          }
3091                  }
3092                  no_glaresources=i;
3093 <
3093 > /*printf("%i",no_glaresources );*/
3094   /* glare sources */
2686        if (detail_out == 1) {
3095  
3096 +        if (output == 4 ) {
3097 +
3098 +        i=0;
3099 +        for (x = 0; x <= igs; x++) {
3100 +                                if (pict_get_npix(p, x) > 0) {
3101 +                                        i = i + 1;
3102 +
3103 +                                        x2=pict_get_av_posx(p, x);
3104 +                                        y2=pict_get_av_posy(p, x);
3105 +
3106 +                                        printf("%i %f %f %f %f %.10f %f %f %f %f \n",
3107 +                                                   i, pict_get_npix(p, x), pict_get_av_posx(p, x), pict_get_av_posy(p, x),
3108 +                                                   pict_get_av_lum(p, x), pict_get_av_omega(p, x),
3109 +                                                   get_posindex(p, pict_get_av_posx(p, x),pict_get_av_posy(p, x),posindex_2),
3110 +                                                   pict_get_cached_dir(p, x2,y2)[0],pict_get_cached_dir(p, x2,y2)[1],pict_get_cached_dir(p, x2,y2)[2]);
3111 +        }
3112 +        }
3113 +
3114 +
3115 +                for (y = 0; y < pict_get_ysize(p); y++) {
3116 +                        for (x = 0; x < pict_get_xsize(p); x++) {
3117 +
3118 +                        if (pict_get_hangle(p, x, y, p->view.vdir, p->view.vup, &ang)) {
3119 +                                if (pict_is_validpixel(p, x, y)) {
3120 +                        if (pict_get_gsn(p,x,y) >0 ) {
3121 +                        i=pict_get_gsn(p,x,y);
3122 +                        printf("%i %i %f %f %f \n", i, new_gs_number[i], pict_get_cached_dir(p, x,y)[0],pict_get_cached_dir(p, x,y)[1],pict_get_cached_dir(p, x,y)[2] );
3123 +                        
3124 +                        }
3125 +
3126 +
3127 +
3128 + }}}}
3129 +
3130 +
3131 + }
3132 +
3133 +
3134 +
3135 +
3136 +        if (detail_out == 1 && output < 3) {
3137 +        dgp = get_dgp(p, E_v, igs, a1, a2, a3, a4, a5, c1, c2, c3, posindex_2);
3138 +
3139                  printf
3140 <                        ("%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",
3140 >                        ("%i No pixels x-pos y-pos L_s Omega_s Posindx L_b L_t E_v Edir Max_Lum Sigma xdir ydir zdir Eglare Lveil_cie teta glare_zone r_contrast r_dgp\n",
3141                           i);
3142                  if (i == 0) {
3143 <                        printf("%i %f %f %f %f %.10f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f\n", i, 0.0, 0.0,
3144 <                                   0.0, 0.0, 0.0, 0.0, lum_backg, lum_task, E_v, E_v_dir,abs_max,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0);
3143 >                        printf("%i %f %f %f %f %.10f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f\n", i, 0.0, 0.0,
3144 >                                   0.0, 0.0, 0.0, 0.0, lum_backg, lum_task, E_v, E_v_dir,abs_max,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0);
3145  
3146                  } else {
3147                          i = 0;
# Line 2710 | Line 3161 | if (calcfast ==1 || search_pix <= 1.0 || calcfast == 2
3161                                                                       Lveil_cie =0 ;
3162                                                                     }
3163                                          Lveil_cie_sum =  Lveil_cie_sum + Lveil_cie;
3164 +                                        r_glare= c2*log10(1 + pow(pict_get_av_lum(p, x),2)* pict_get_av_omega(p, x)/pow(get_posindex(p, pict_get_av_posx(p, x), pict_get_av_posy(p, x),posindex_2), 2.0) / pow(E_v, a3))  ;
3165 +                                        r_contrast=r_glare/sum_glare  ;
3166 +                                        r_dgp=r_glare/dgp ;
3167                                          if (pict_get_z1_gsn(p,x)<0) {
3168                                          act_gsn=(int)(-pict_get_z1_gsn(p,x));
3169                                          }else{
3170                                          act_gsn=0;
3171                                          }
3172 <                                        printf("%i %f %f %f %f %.10f %f %f %f %f %f %f %f %f %f %f %f %f %f %i \n",
3172 >                                        printf("%i %f %f %f %f %.10f %f %f %f %f %f %f %f %f %f %f %f %f %f %i %f %f \n",
3173                                                     i, pict_get_npix(p, x), pict_get_av_posx(p, x),
3174                                                     pict_get_ysize(p) - pict_get_av_posy(p, x),
3175                                                     pict_get_av_lum(p, x), pict_get_av_omega(p, x),
# Line 2723 | Line 3177 | if (calcfast ==1 || search_pix <= 1.0 || calcfast == 2
3177                                                                                  pict_get_av_posy(p, x),
3178                                                                                  posindex_2), lum_backg, lum_task,
3179                                                     E_v, E_v_dir, abs_max, sigma * 180 / 3.1415927
3180 <                                                   ,pict_get_cached_dir(p, x2,y2)[0],pict_get_cached_dir(p, x2,y2)[1],pict_get_cached_dir(p, x2,y2)[2],pict_get_Eglare(p,x),Lveil_cie,teta,act_gsn );
3180 >                                                   ,pict_get_cached_dir(p, x2,y2)[0],pict_get_cached_dir(p, x2,y2)[1],pict_get_cached_dir(p, x2,y2)[2],pict_get_Eglare(p,x),Lveil_cie,teta,act_gsn,r_contrast,r_dgp );
3181                                  }
3182                          }
3183                  }
3184          }
3185  
3186 + if ( output < 3) {
3187  
2733
3188   /* calculation of indicees */
3189  
3190   /* check vertical illuminance range */
# Line 2743 | Line 3197 | if (calcfast ==1 || search_pix <= 1.0 || calcfast == 2
3197          dgp =
3198                  get_dgp(p, E_v2, igs, a1, a2, a3, a4, a5, c1, c2, c3, posindex_2);
3199   /* low light correction */
3200 +     if (lowlight ==1) {
3201         if (E_v < 1000) {
3202         low_light_corr=1.0*exp(0.024*E_v-4)/(1+exp(0.024*E_v-4));} else {low_light_corr=1.0 ;}
3203         dgp =low_light_corr*dgp;
3204 <      
3204 >       }
3205   /* age correction */
3206        
3207          if (age_corr == 1) {
# Line 2767 | Line 3222 | if (calcfast == 0) {
3222          dgi = get_dgi(p, lum_backg, igs, posindex_2);
3223          ugr = get_ugr(p, lum_backg, igs, posindex_2);
3224          ugp = get_ugp (p,lum_backg, igs, posindex_2);
3225 +        ugp2 = get_ugp2 (p,lum_backg, igs, posindex_2);
3226          ugr_exp = get_ugr_exp (p,lum_backg_cos,lum_a, igs, posindex_2);
3227          dgi_mod = get_dgi_mod(p, lum_a, igs, posindex_2);
3228          cgi = get_cgi(p, E_v, E_v_dir, igs, posindex_2);
# Line 2777 | Line 3233 | if (calcfast == 0) {
3233                  dgi = 0.0;
3234                  ugr = 0.0;
3235                  ugp = 0.0;
3236 +                ugp2 =0.0;
3237                  ugr_exp = 0.0;
3238                  dgi_mod = 0.0;
3239                  cgi = 0.0;
# Line 2803 | Line 3260 | if (calcfast == 0) {
3260                                 get_dgp(p, E_vl_ext, igs, a1, a2, a3, a4, a5, c1, c2, c3,
3261                                  posindex_2);
3262                                  dgp = dgp_ext;
3263 <                                if (E_vl_ext < 1000) {
3263 >                                if (E_vl_ext < 1000 && lowlight ==1) {
3264                                  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 ;}
3265                                  dgp =low_light_corr*dgp;
3266                                  dgp =age_corr_factor*dgp;
# Line 2816 | Line 3273 | if (calcfast == 0) {
3273  
3274          
3275                          printf
3276 <                                ("dgp,av_lum,E_v,lum_backg,E_v_dir,dgi,ugr,vcp,cgi,lum_sources,omega_sources,Lveil,Lveil_cie,dgr,ugp,ugr_exp,dgi_mod,av_lum_pos,av_lum_pos2,med_lum,med_lum_pos,med_lum_pos2: %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f\n",
3276 >                                ("dgp,av_lum,E_v,lum_backg,E_v_dir,dgi,ugr,vcp,cgi,lum_sources,omega_sources,Lveil,Lveil_cie,dgr,ugp,ugr_exp,dgi_mod,av_lum_pos,av_lum_pos2,med_lum,med_lum_pos,med_lum_pos2,ugp2: %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f\n",
3277                                   dgp, avlum, E_v, lum_backg, E_v_dir, dgi, ugr, vcp, cgi,
3278 <                                 lum_sources, omega_sources, Lveil,Lveil_cie_sum,dgr,ugp,ugr_exp,dgi_mod,lum_pos_mean,lum_pos2_mean/sang,lum_nopos_median[0],lum_pos_median[0],lum_pos2_median[0]);
3278 >                                 lum_sources, omega_sources, Lveil,Lveil_cie_sum,dgr,ugp,ugr_exp,dgi_mod,lum_pos_mean,lum_pos2_mean/sang,lum_nopos_median[0],lum_pos_median[0],lum_pos2_median[0],ugp2);
3279                  } else {
3280                          if (detail_out2 == 1) {
3281  
# Line 2831 | Line 3288 | if (calcfast == 0) {
3288                                  if (ext_vill == 1) {
3289                  dgp_ext = get_dgp(p, E_vl_ext, igs, a1, a2, a3, a4, a5, c1, c2, c3,posindex_2);
3290                                  
3291 <                                if (E_vl_ext < 1000) {
3291 >                                if (E_vl_ext < 1000 && lowlight ==1) {
3292                                  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 ;}
3293                                          dgp =low_light_corr*dgp_ext;
3294                                          dgp =age_corr_factor*dgp;
# Line 2846 | Line 3303 | if (calcfast == 0) {
3303                                 get_dgp(p, E_vl_ext, igs, a1, a2, a3, a4, a5, c1, c2, c3,
3304                                  posindex_2);
3305                                  dgp = dgp_ext;
3306 <                                if (E_vl_ext < 1000) {
3306 >                                if (E_vl_ext < 1000 && lowlight ==1) {
3307                                  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 ;}
3308                                  dgp =low_light_corr*dgp;
3309  
# Line 2859 | Line 3316 | if (calcfast == 0) {
3316                           printf("%f\n", dgp);
3317                  }
3318          }
3319 + }
3320  
3321 + }else{
3322 + /* only multiimagemode                        
3323 +
3324 + */
3325 +
3326 +
3327 +                       for (j = 0; j < num_images; j++) {
3328 +
3329 +                                
3330 + /* loop over temporal images */
3331 +
3332 + pict_read(pcoeff,temp_image_name[j] );
3333 +
3334 + /*printf ("num_img:%i x-size:%i xpos:%i y-size: %i ypos %i bigimg-xsize %i %i ",num_images,pict_get_xsize(pcoeff),x_temp_img[j],pict_get_ysize(pcoeff),y_temp_img[j],pict_get_xsize(p),pict_get_ysize(p));
3335 + */
3336 +
3337 +
3338 + /* loop over scans
3339 +    empty of */
3340 + for (jj=0;jj<=num_scans;jj++){
3341 +
3342 + /*printf("scale %i %i %i %f ",j,jj,num_scans,scale_image_scans[j*(num_scans+1)+jj]); */
3343 +
3344 +
3345 + /* copy luminance value into big image and remove glare sources*/
3346 +        for (x = 0; x < pict_get_xsize(pcoeff); x++) {
3347 +                for (y = 0; y < pict_get_ysize(pcoeff); y++) {
3348 +                        xmap=x_temp_img[j]+x;
3349 +                        ymap=y_temp_img[j]+y;
3350 +                        if (xmap <0) { xmap=0;}
3351 +                        if (ymap <0) { ymap=0;}
3352 +                        
3353 +                        pict_get_color(p, xmap, ymap)[RED] = scale_image_scans[j*(num_scans+1)+jj]*pict_get_color(pcoeff, x, y)[RED];
3354 +                        pict_get_color(p, xmap, ymap)[GRN] = scale_image_scans[j*(num_scans+1)+jj]*pict_get_color(pcoeff, x, y)[GRN];
3355 +                        pict_get_color(p, xmap, ymap)[BLU] = scale_image_scans[j*(num_scans+1)+jj]*pict_get_color(pcoeff, x, y)[BLU];
3356 +                        pict_get_gsn(p, xmap, ymap) = 0;
3357 +                        pict_get_pgs(p, xmap, ymap) = 0;
3358 + }}
3359 +
3360 +
3361 +
3362 +
3363 +
3364 + actual_igs =0;
3365 +
3366 + /* first glare source scan: find primary glare sources */
3367 +        for (x = 0; x < pict_get_xsize(pcoeff); x++) {
3368 +                lastpixelwas_gs=0;
3369 + /*              lastpixelwas_peak=0; */
3370 +                for (y = 0; y < pict_get_ysize(pcoeff); y++) {
3371 +                        xmap=x_temp_img[j]+x;
3372 +                        ymap=y_temp_img[j]+y;
3373 +                        if (pict_get_hangle(p, xmap, ymap, p->view.vdir, p->view.vup, &ang)) {
3374 +                                if (pict_is_validpixel(p, xmap, ymap)) {
3375 +                                        act_lum = luminance(pict_get_color(p, xmap, ymap));
3376 +                                        if (act_lum> lum_source) {
3377 +                                                if (act_lum >lum_total_max) {
3378 +                                                                             lum_total_max=act_lum;
3379 +                                                                               }
3380 +                      
3381 +                                                if (lastpixelwas_gs==0 || search_pix <= 1.0 ) {
3382 +                                                actual_igs = find_near_pgs(p, xmap, ymap, max_angle, 0, igs, 1);
3383 +  }
3384 +                                                if (actual_igs == 0) {
3385 +                                                        igs = igs + 1;
3386 +                                                        pict_new_gli(p);
3387 +                                                        pict_get_Eglare(p,igs) = 0.0;
3388 + /*  not necessary here                                  pict_get_lum_min(p, igs) = HUGE_VAL;
3389 +                                                        pict_get_Eglare(p,igs) = 0.0;
3390 +                                                        pict_get_Dglare(p,igs) = 0.0;
3391 +                                                        pict_get_z1_gsn(p,igs) = 0;
3392 +                                                        pict_get_z2_gsn(p,igs) = 0; */
3393 +                                                        actual_igs = igs;
3394 +                                                        
3395 +                                                }
3396 +                                                pict_get_gsn(p, xmap, ymap) = actual_igs;
3397 +                                                pict_get_pgs(p, xmap, ymap) = 1;
3398 +                                                add_pixel_to_gs(p, xmap, ymap, actual_igs);
3399 +                                                lastpixelwas_gs=actual_igs;
3400 +                                                
3401 +                                                
3402 +                                                
3403 +                                                delta_E=DOT(p->view.vdir, pict_get_cached_dir(p, xmap, ymap))* pict_get_omega(p, xmap, ymap)* luminance(pict_get_color(p, xmap, ymap));
3404 +                                                pict_get_Eglare(p,actual_igs ) = pict_get_Eglare(p,actual_igs ) + delta_E;
3405 +
3406 +                                              
3407 +                                                
3408 +
3409 +                                        } else {
3410 +                                                pict_get_pgs(p, xmap, ymap) = 0;
3411 +                                                lastpixelwas_gs=0;
3412 +                                        }
3413 +                                }
3414 +                        }
3415 +                }
3416 +             }
3417 +
3418 +
3419 + /*                              here should be peak extraction  */
3420 + i_max=igs;
3421 +                r_split = max_angle / 2.0;
3422 +                for (i = 0; i <= i_max; i++) {
3423 +
3424 +                        if (pict_get_npix(p, i) > 0) {
3425 +                                l_max = pict_get_lum_max(p, i);
3426 +                                i_splitstart = igs + 1;
3427 +                                if (l_max >= limit) {
3428 +                                        for (x = 0; x < pict_get_xsize(pcoeff); x++)
3429 +                                                for (y = 0; y < pict_get_ysize(pcoeff); y++) {
3430 +                                                xmap=x_temp_img[j]+x;
3431 +                                                ymap=y_temp_img[j]+y;
3432 +                                                
3433 +                                                
3434 +                                                        if (pict_get_hangle
3435 +                                                                (p, xmap, ymap, p->view.vdir, p->view.vup, &ang))
3436 +                                                        {
3437 +                                                                if (pict_is_validpixel(p, xmap, ymap)
3438 +                                                                        && luminance(pict_get_color(p, xmap, ymap))
3439 +                                                                        >= limit
3440 +                                                                        && pict_get_gsn(p, xmap, ymap) == i) {
3441 +                                                                        if (i_splitstart == (igs + 1)) {
3442 +                                                                                igs = igs + 1;
3443 +                                                                                pict_new_gli(p);
3444 +                                                                                pict_get_z1_gsn(p,igs) = 0;
3445 +                                                                                pict_get_z2_gsn(p,igs) = 0;
3446 +
3447 +                                                                                pict_get_Eglare(p,igs) = 0.0;
3448 +                                                                                pict_get_Dglare(p,igs) = 0.0;
3449 +                                                                                pict_get_lum_min(p, igs) =
3450 +                                                                                        99999999999999.999;
3451 +                                                                                i_split = igs;
3452 +                                                                        } else {
3453 +                                                                                i_split =
3454 +                                                                                        find_split(p, xmap, ymap, r_split,
3455 +                                                                                                           i_splitstart, igs);
3456 +                                                                        }
3457 +                                                                        if (i_split == 0) {
3458 +                                                                                igs = igs + 1;
3459 +                                                                                pict_new_gli(p);
3460 +                                                                                pict_get_z1_gsn(p,igs) = 0;
3461 +                                                                                pict_get_z2_gsn(p,igs) = 0;
3462 +
3463 +                                                                                pict_get_Eglare(p,igs) = 0.0;
3464 +                                                                                pict_get_Dglare(p,igs) = 0.0;
3465 +                                                                                pict_get_lum_min(p, igs) =
3466 +                                                                                        99999999999999.999;
3467 +                                                                                i_split = igs;
3468 +                                                                        }
3469 +                                                                        split_pixel_from_gs(p, xmap, ymap, i_split, uniform_gs, u_r, u_g , u_b);
3470 +
3471 +                                                                }
3472 +                                                        }
3473 +                                                }
3474 +
3475 +                                }
3476 +                                change = 1;
3477 +                                while (change == 1) {
3478 +                                        change = 0;
3479 +                                        for (x = 0; x < pict_get_xsize(pcoeff); x++)
3480 +                                                for (y = 0; y < pict_get_ysize(pcoeff); y++) {
3481 +                                                xmap=x_temp_img[j]+x;
3482 +                                                ymap=y_temp_img[j]+y;
3483 +                                                        before_igs = pict_get_gsn(p, xmap, ymap);
3484 +                                                        if (before_igs >= i_splitstart) {
3485 +                                                                if (pict_get_hangle
3486 +                                                                        (p, xmap, ymap, p->view.vdir, p->view.vup,
3487 +                                                                         &ang)) {
3488 +                                                                        if (pict_is_validpixel(p, xmap, ymap)
3489 +                                                                                && before_igs > 0) {
3490 +                                                                                actual_igs =
3491 +                                                                                        find_near_pgs(p, xmap, ymap,
3492 +                                                                                                                  max_angle,
3493 +                                                                                                                  before_igs, igs,
3494 +                                                                                                                  i_splitstart);
3495 +                                                                                if (!(actual_igs == before_igs)) {
3496 +                                                                                        change = 1;
3497 +                                                                                }
3498 +                                                                                if (before_igs > 0) {
3499 +                                                                                        actual_igs =
3500 +                                                                                                pict_get_gsn(p, xmap, ymap);
3501 + /*                                                                                     icol =
3502 +                                                                                                setglcolor(p, x, y,
3503 +                                                                                                                   actual_igs, uniform_gs, u_r, u_g , u_b);*/
3504 +                                                                                }
3505 +
3506 +                                                                        }
3507 +                                                                }
3508 +                                                        }
3509 +
3510 +                                                }
3511 +                                }
3512 +
3513 +
3514 +                        }
3515 +                }
3516 +
3517 + /*                              end peak extraction  */
3518 +
3519 +
3520 + /* calculation of direct vertical illuminance for th multi-image-mode */
3521 +
3522 +
3523 +        for (x = 0; x < pict_get_xsize(pcoeff); x++)
3524 +                for (y = 0; y < pict_get_ysize(pcoeff); y++) {
3525 +                        xmap=x_temp_img[j]+x;
3526 +                        ymap=y_temp_img[j]+y;
3527 +                        if (pict_get_hangle(p, xmap, ymap, p->view.vdir, p->view.vup, &ang)) {
3528 +                                if (pict_is_validpixel(p, xmap, ymap)) {
3529 +                                        if (pict_get_gsn(p, xmap, ymap) > 0) {
3530 +                                                actual_igs = pict_get_gsn(p, xmap, ymap);
3531 +                                                delta_E=DOT(p->view.vdir, pict_get_cached_dir(p, xmap, ymap))* pict_get_omega(p, xmap, ymap)* luminance(pict_get_color(p, xmap, ymap));
3532 +                                                pict_get_Eglare(p,actual_igs ) = pict_get_Eglare(p,actual_igs ) + delta_E;
3533 +                                        }
3534 +                                }
3535 +                        }
3536 +                }
3537 +
3538 +
3539 +
3540 +
3541 +
3542 +
3543 +
3544 +
3545 +                        i = 0;
3546 +                        for (x = 0; x <= igs; x++) {
3547 +                                if (pict_get_npix(p, x) > 0) {
3548 +                                        i = i + 1;
3549 +                                        }}
3550 + no_glaresources=i;                      
3551 +
3552 + /*
3553 +
3554 + sprintf(file_out, "%s%i%s","ray2/img_",j,".hdr");
3555 + pict_write(p, file_out);
3556 + */
3557 + printf("%i ",no_glaresources);
3558 +                        i = 0;
3559 +                        for (x = 0; x <= igs; x++) {
3560 +                                if (pict_get_npix(p, x) > 0) {
3561 +                                        i = i + 1;
3562 +                                        pict_get_sigma(p, pict_get_av_posx(p, x), pict_get_av_posy(p, x), p->view.vdir, p->view.vup, &sigma);
3563 +                                                                  
3564 +                                        x2=pict_get_av_posx(p, x);
3565 +                                        y2=pict_get_av_posy(p, x);
3566 +                                        printf("%f %.10f %f %f %f %f %f ",      pict_get_av_lum(p, x), pict_get_av_omega(p, x), get_posindex(p, pict_get_av_posx(p, x), pict_get_av_posy(p, x),
3567 +                                                                posindex_2),pict_get_cached_dir(p, x2,y2)[0],pict_get_cached_dir(p, x2,y2)[1],pict_get_cached_dir(p, x2,y2)[2],pict_get_Eglare(p,x) );
3568 +                                }
3569 + pict_get_npix(p, x)=0;
3570 + pict_get_av_lum(p, x)=0;
3571 + pict_get_av_posy(p, x)=0;
3572 + pict_get_av_posx(p, x)=0;
3573 + pict_get_av_omega(p, x)=0;
3574 +                        }
3575 + printf("\n");
3576 +
3577 +
3578 + /* empty big image and remove glare sources*/
3579 +        for (x = 0; x < pict_get_xsize(pcoeff); x++) {
3580 +                for (y = 0; y < pict_get_ysize(pcoeff); y++) {
3581 +                        xmap=x_temp_img[j]+x;
3582 +                        ymap=y_temp_img[j]+y;
3583 +                        pict_get_color(p, xmap, ymap)[RED] = 0;
3584 +                        pict_get_color(p, xmap, ymap)[GRN] = 0;
3585 +                        pict_get_color(p, xmap, ymap)[BLU] = 0;
3586 +                        pict_get_gsn(p, xmap, ymap) = 0;
3587 +                        pict_get_pgs(p, xmap, ymap) = 0;
3588 + }}
3589 + igs=0;
3590 +
3591 +
3592 + }
3593 +
3594 +
3595 + }
3596 +
3597 + }
3598 +
3599 + /* end multi-image-mode */
3600  
3601   /*printf("lowlight=%f\n", low_light_corr); */
3602  

Diff Legend

Removed lines
+ Added lines
< Changed lines (old)
> Changed lines (new)