ViewVC Help

View File | Revision Log | Show Annotations | Download File | Root Listing

root/radiance/ray/src/util/evalglare.c

Comparing ray/src/util/evalglare.c (file contents):
Revision 2.1 by greg, Wed Aug 12 23:07:59 2015 UTC vs.
Revision 2.9 by greg, Fri Nov 23 19:32:17 2018 UTC

#	Line 1 | Line 1
1	<	/* EVALGLARE V1.17
2	<	* Evalglare Software License, Version 1.0
1	>	#ifndef lint
2	>	static const char RCSid[] = "$Id$";
3	>	#endif
4	>	/* EVALGLARE V2.07
5	>	* Evalglare Software License, Version 2.0
6		*
7	<	* Copyright (c) 1995 - 2015 Fraunhofer ISE, EPFL.
7	>	* Copyright (c) 1995 - 2016 Fraunhofer ISE, EPFL.
8		* All rights reserved.
9		*
10		*
#	Line 20 | Line 23
23		* 3. The end-user documentation included with the redistribution,
24		*           if any, must include the following acknowledgments:
25		*             "This product includes the evalglare software,
26	<	developed at Fraunhofer ISE by J. Wienold" and
26	>	developed at Fraunhofer ISE and EPFL by J. Wienold" and
27		*             "This product includes Radiance software
28		*                 (http://radsite.lbl.gov/)
29		*                 developed by the Lawrence Berkeley National Laboratory
#	Line 28 | Line 31
31		*       Alternately, this acknowledgment may appear in the software itself,
32		*       if and wherever such third-party acknowledgments normally appear.
33		*
34	<	* 4. The names "Evalglare," and "Fraunhofer ISE" must
34	>	* 4. The names "Evalglare," "EPFL" and "Fraunhofer ISE" must
35		*       not be used to endorse or promote products derived from this
36		*       software without prior written permission. For written
37		*       permission, please contact [email protected]
38		*
39		* 5. Products derived from this software may not be called "evalglare",
40		*       nor may "evalglare" appear in their name, without prior written
41	<	*       permission of Fraunhofer ISE.
41	>	*       permission of EPFL.
42		*
43		* Redistribution and use in source and binary forms, WITH
44		* modification, are permitted provided that the following conditions
#	Line 45 | Line 48
48		* conditions 1.-5. apply
49		*
50		* 6. In order to ensure scientific correctness, any modification in source code imply fulfilling all following comditions:
51	<	*     a) A written permission from Fraunhofer ISE. For written permission, please contact [email protected].
51	>	*     a) A written permission from EPFL. For written permission, please contact [email protected].
52		*     b) The permission can be applied via email and must contain the applied changes in source code and at least two example calculations,
53		*        comparing the results of the modified version and the current version of evalglare.
54		*     b) Any redistribution of a modified version of evalglare must contain following note:
#	Line 275 | Line 278
278		remove of age factor due to non proven statistical evidence
279		*/
280
281	<	#define EVALGLARE
281	>	/* evalglare.c, v1.18 2015/11/10 add masking, UDP, PSGV, DGI_mod, UGR_exp and zoning. For zoning and band luminance, a statistical package from C. Reetz was added in order to derive median, std_dev, percentiles.
282	>	memory leakage was checked and is o.k.
283	>	*/
284
285	+	/* evalglare.c, v1.19 2015/12/04 sorting subroutines in statistical package changed by. C. Reetz, adding Position index weighted Luminance evaluation (mean and median of entire image, only in detailed output available)
286	+	*/
287	+	/* evalglare.c, v1.20 2015/12/21 removing -ansi gcc-option in makefile for the standalone version, adding DGR as output since VCP equation seems not to match with original publication
288	+	*/
289	+	/* evalglare.c, v1.21 2016/03/15 add a new pixel-overflow-correction option -N
290	+	*/
291	+	/* evalglare.c, v1.22 2016/03/22 fixed problem using -N option and angular distance calculation, when actual pixel equals disk-center (arccos returned nan)
292	+	*/
293	+	/* evalglare.c, v1.23 2016/04/18 fixed problem making black corners for -vth fish-eye images
294	+	*/
295	+	/* evalglare.c, v1.24 2016/04/26 significant speed improvement for 2nd glare source scan: for other glare sources will not be searched any more, when the 8 surrounding pixels have the same glare source number! Smaller speed improvement in the first glare source scan: remembering if pixel before was also a glare source, then no search for surrounding glare sources is applied.
296	+	changed masking threshold to 0.05 cd/m2
297	+	*/
298	+	/* evalglare.c, v1.25 2016/04/27 removed the change on the first glare source scan: in few cases glare sources are not merged together in the same way as before. Although the difference is marginal, this algorithm was removed in order to be consistent to the existing results.
299	+	*/
300	+	/* evalglare.c, v1.26 2016/04/28 removed the bug Lb getting nan in case all glare source pixels are above the peak extraction limit.
301	+	The calculation of the average source lumiance and solid angle was adapted to be more robust for extreme cases.
302	+	In case no glare source is found, the values of the glare metrics, relying only on glare sources, is set to zero (100 for VCP) to avoid nan and -inf in the output.
303	+	Changed glare section output when 0 glare source are found to have the same amount of columns than when glare sources are found
304	+	*/
305	+	/* evalglare.c, v1.27 2016/06/13 include a warning, if -vtv is in the header. Check, if the corners are black AND -vtv is used ->stopping (stopping can be avoided by using the forcing option -f ). Check, if an invalid exposure is in the header -> stopping in any case.
306	+	*/
307	+	/* evalglare.c, v1.28 2016/07/08 small code changes in the section of calculating the E_glare (disability) for each glare source.
308	+	*/
309	+	/* evalglare.c, v1.29 2016/07/12 change threshold for the black corner to 2cd/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.17 release 15.07.2015 by EPFL, J.Wienold"
346	>	#define VERSION "2.07 release 17.11.2018 by EPFL, J.Wienold"
347		#define RELEASENAME PROGNAME " " VERSION
348
349	<	#include "rtio.h"
285	<	#include "platform.h"
349	>
350		#include "pictool.h"
351	+	#include "rtio.h"
352		#include <math.h>
353		#include <string.h>
354	+	#include "platform.h"
355	+	#include "muc_randvar.h"
356	+
357		char *progname;
358
359		/* subroutine to add a pixel to a glare source */
#	Line 561 | Line 629 | double get_task_lum(pict * p, int x, int y, float r, i
629		return av_lum;
630		}
631
564	–	/* subroutine for calculation of band luminance */
565	–	double get_band_lum(pict * p, float r, int task_color)
566	–	{
567	–	int x_min, x_max, y_min, y_max, ix, iy, y_mid;
568	–	double r_actual, av_lum, omega_sum, act_lum;
632
633
571	–	x_max = pict_get_xsize(p) - 1;
572	–	y_max = pict_get_ysize(p) - 1;
573	–	x_min = 0;
574	–	y_min = 0;
575	–	y_mid = rint(y_max/2);
634
635
578	–
579	–	av_lum = 0.0;
580	–	omega_sum = 0.0;
581	–
582	–	for (iy = y_min; iy <= y_max; iy++) {
583	–	for (ix = x_min; ix <= x_max; ix++) {
584	–
585	–	/*                      if (DOT(pict_get_cached_dir(p,ix,iy),p->view.vdir) < 0.0)
586	–	continue;*/
587	–	r_actual =
588	–	acos(DOT(pict_get_cached_dir(p, ix, y_mid),
589	–	pict_get_cached_dir(p, ix, iy))) ;
590	–	act_lum = luminance(pict_get_color(p, ix, iy));
591	–
592	–	if ((r_actual <= r) || (iy == y_mid) ) {
593	–	act_lum = luminance(pict_get_color(p, ix, iy));
594	–	av_lum += pict_get_omega(p, ix, iy) * act_lum;
595	–	omega_sum += pict_get_omega(p, ix, iy);
596	–	if (task_color == 1) {
597	–	pict_get_color(p, ix, iy)[RED] = 0.0;
598	–	pict_get_color(p, ix, iy)[GRN] = act_lum / WHTEFFICACY / CIE_gf;
599	–	pict_get_color(p, ix, iy)[BLU] = 0.0;
600	–	}
601	–	}
602	–	}
603	–	}
604	–
605	–	av_lum = av_lum / omega_sum;
606	–	/*    printf("average luminance of band %f \n",av_lum);*/
607	–	/*      printf("solid angle of band %f \n",omega_sum);*/
608	–	return av_lum;
609	–	}
610	–
611	–
612	–
613	–
614	–
636		/* subroutine for coloring the glare sources
637		red is used only for explicit call of the subroutine with acol=0  , e.g. for disability glare
638		-1: set to grey*/
#	Line 623 | Line 644 | int setglcolor(pict * p, int x, int y, int acol, int u
644		l=u_r+u_g+u_b ;
645
646		pcol[RED][0] = 1.0 / CIE_rf;
647	<	pcol[GRN][0] = 0.;
648	<	pcol[BLU][0] = 0.;
647	>	pcol[GRN][0] = 0.0 / CIE_gf;
648	>	pcol[BLU][0] = 0.0 / CIE_bf;
649
650	<	pcol[RED][1] = 0.;
650	>	pcol[RED][1] = 0.0 / CIE_rf;
651		pcol[GRN][1] = 1.0 / CIE_gf;
652	<	pcol[BLU][1] = 0.;
652	>	pcol[BLU][1] = 0.0 / CIE_bf;
653
654	<	pcol[RED][2] = 0.;
655	<	pcol[GRN][2] = 0.;
656	<	pcol[BLU][2] = 1 / CIE_bf;
654	>	pcol[RED][2] = 0.15 / CIE_rf;
655	>	pcol[GRN][2] = 0.15 / CIE_gf;
656	>	pcol[BLU][2] = 0.7 / CIE_bf;
657
658	<	pcol[RED][3] = 1.0 / (1.0 - CIE_bf);
659	<	pcol[GRN][3] = 1.0 / (1.0 - CIE_bf);
660	<	pcol[BLU][3] = 0.;
658	>	pcol[RED][3] = .5 / CIE_rf;
659	>	pcol[GRN][3] = .5 / CIE_gf;
660	>	pcol[BLU][3] = 0.0 / CIE_bf;
661
662	<	pcol[RED][4] = 1.0 / (1.0 - CIE_gf);
663	<	pcol[GRN][4] = 0.;
664	<	pcol[BLU][4] = 1.0 / (1.0 - CIE_gf);
662	>	pcol[RED][4] = .5 / CIE_rf;
663	>	pcol[GRN][4] = .0 / CIE_gf;
664	>	pcol[BLU][4] = .5 / CIE_bf ;
665
666	<	pcol[RED][5] = 0.;
667	<	pcol[GRN][5] = 1.0 / (1.0 - CIE_rf);
668	<	pcol[BLU][5] = 1.0 / (1.0 - CIE_rf);
666	>	pcol[RED][5] = .0 / CIE_rf;
667	>	pcol[GRN][5] = .5 / CIE_gf;
668	>	pcol[BLU][5] = .5 / CIE_bf;
669
670	<	pcol[RED][6] = 0.;
671	<	pcol[GRN][6] = 1.0 / (1.0 - CIE_rf);
672	<	pcol[BLU][6] = 1.0 / (1.0 - CIE_rf);
670	>	pcol[RED][6] = 0.333 / CIE_rf;
671	>	pcol[GRN][6] = 0.333 / CIE_gf;
672	>	pcol[BLU][6] = 0.333 / CIE_bf;
673
674
675		pcol[RED][999] = 1.0 / WHTEFFICACY;
#	Line 659 | Line 680 | int setglcolor(pict * p, int x, int y, int acol, int u
680		pcol[RED][998] = u_r /(l* CIE_rf) ;
681		pcol[GRN][998] = u_g /(l* CIE_gf);
682		pcol[BLU][998] = u_b /(l* CIE_bf);
683	<	/*        printf("CIE_rf,CIE_gf,CIE_bf,l=%f %f %f %f\n",CIE_rf,CIE_gf,CIE_bf,l);*/
683	>	/*        printf("CIE_rf,CIE_gf,CIE_bf,l=%f %f %f %f %f\n",CIE_rf,CIE_gf,CIE_bf,l,WHTEFFICACY); */
684		icol = acol ;
685		if ( acol == -1) {icol=999;
686		}else{if (acol>0){icol = acol % 5 +1;
#	Line 672 | Line 693 | int setglcolor(pict * p, int x, int y, int acol, int u
693		pict_get_color(p, x, y)[RED] = pcol[RED][icol] * act_lum / WHTEFFICACY;
694		pict_get_color(p, x, y)[GRN] = pcol[GRN][icol] * act_lum / WHTEFFICACY;
695		pict_get_color(p, x, y)[BLU] = pcol[BLU][icol] * act_lum / WHTEFFICACY;
696	<
696	>	/*        printf("x,y,lum_before,lum_after,diff %i %i %f %f %f\n",x,y, act_lum,luminance(pict_get_color(p, x, y)), act_lum-luminance(pict_get_color(p, x, y))); */
697		return icol;
698		}
699
#	Line 724 | Line 745 | void add_secondary_gs(pict * p, int x, int y, float r,
745
746		/* color pixel of gs */
747
748	<	icol = setglcolor(p, x, y, act_gs, uniform_gs, u_r, u_g , u_b);
748	>	/*              icol = setglcolor(p, x, y, act_gs, uniform_gs, u_r, u_g , u_b); */
749
750		}
751		}
#	Line 764 | Line 785 | void split_pixel_from_gs(pict * p, int x, int y, int n
785
786		/* color pixel of new gs */
787
788	<	icol = setglcolor(p, x, y, new_gsn, uniform_gs, u_r, u_g , u_b);
788	>	/*      icol = setglcolor(p, x, y, new_gsn, uniform_gs, u_r, u_g , u_b); */
789	>
790
769	–
791		}
792
793		/* subroutine for the calculation of the position index */
#	Line 797 | Line 818 | float get_posindex(pict * p, float x, float y, int pos
818		tau = tau * deg;
819		sigma = sigma * deg;
820
821	+	if (postype == 1) {
822	+	/* KIM  model */
823	+	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));
824	+	}else{
825	+	/* Guth model, equation from IES lighting handbook */
826		posindex =
827		exp((35.2 - 0.31889 * tau -
828		1.22 * exp(-2 * tau / 9)) / 1000 * sigma + (21 +
#	Line 818 | Line 844 | float get_posindex(pict * p, float x, float y, int pos
844
845		posindex = 1 + fact * r;
846		}
847	<	if (posindex > 16)
847	>	if (posindex > 16)
848		posindex = 16;
849	+	}
850
851		return posindex;
852
#	Line 1032 | Line 1059 | return;
1059
1060		}
1061
1062	<	/* subroutine for the calculation of the daylight glare index */
1062	>	/* subroutine for the calculation of the daylight glare index DGI*/
1063
1064
1065		float get_dgi(pict * p, float lum_backg, int igs, int posindex_2)
#	Line 1059 | Line 1086 | float get_dgi(pict * p, float lum_backg, int igs, int
1086		return dgi;
1087
1088		}
1089	+	/* subroutine for the calculation of the modified daylight glare index DGI_mod (Fisekis2003)
1090	+	several glare sources are taken into account and summed up, original equation has no summary
1091	+	*/
1092
1093	+
1094	+	float get_dgi_mod(pict * p, float lum_a, int igs, int posindex_2)
1095	+	{
1096	+	float dgi_mod, sum_glare, omega_s;
1097	+	int i;
1098	+
1099	+
1100	+	sum_glare = 0;
1101	+	omega_s = 0;
1102	+
1103	+	for (i = 0; i <= igs; i++) {
1104	+
1105	+	if (pict_get_npix(p, i) > 0) {
1106	+	omega_s =
1107	+	pict_get_av_omega(p, i) / get_posindex(p,pict_get_av_posx(p,i), pict_get_av_posy(p,i),posindex_2) /
1108	+	get_posindex(p, pict_get_av_posx(p, i), pict_get_av_posy(p, i), posindex_2);
1109	+	sum_glare += 0.478 * pow(pict_get_av_lum(p, i), 1.6) * pow(omega_s,0.8) / (pow(lum_a,0.85) + 0.07 * pow(pict_get_av_omega(p, i),0.5) * pict_get_av_lum(p, i));
1110	+	/*    printf("i,sum_glare %i %f\n",i,sum_glare); */
1111	+	}
1112	+	}
1113	+	dgi_mod = 10 * log10(sum_glare);
1114	+
1115	+	return dgi_mod;
1116	+
1117	+	}
1118	+
1119		/* subroutine for the calculation of the daylight glare probability */
1120		double
1121		get_dgp(pict * p, double E_v, int igs, double a1, double a2, double a3,
#	Line 1100 | Line 1156 | get_dgp(pict * p, double E_v, int igs, double a1, doub
1156		}
1157
1158		/* subroutine for the calculation of the visual comfort probability */
1159	<	float get_vcp(pict * p, double lum_a, int igs, int posindex_2)
1159	>	float get_dgr(pict * p, double lum_a, int igs, int posindex_2)
1160		{
1161	<	float vcp;
1162	<	double sum_glare, dgr;
1161	>	float dgr;
1162	>	double sum_glare;
1163		int i, i_glare;
1164
1165
#	Line 1129 | Line 1185 | float get_vcp(pict * p, double lum_a, int igs, int pos
1185		}
1186		dgr = pow(sum_glare, pow(i_glare, -0.0914));
1187
1188	<	vcp = 50 * erf((6.347 - 1.3227 * log(dgr)) / 1.414213562373) + 50;
1188	>
1189	>
1190	>	return dgr;
1191	>
1192	>	}
1193	>
1194	>	float get_vcp(float dgr )
1195	>	{
1196	>	float vcp;
1197	>
1198	>	vcp = 50 * erf((6.374 - 1.3227 * log(dgr)) / 1.414213562373) + 50;
1199		if (dgr > 750) {
1200		vcp = 0;
1201		}
1202		if (dgr < 20) {
1203		vcp = 100;
1204		}
1139	–
1140	–
1205		return vcp;
1206
1207		}
1208
1209	+
1210		/* subroutine for the calculation of the unified glare rating */
1211		float get_ugr(pict * p, double lum_backg, int igs, int posindex_2)
1212		{
#	Line 1165 | 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		}
1243	+	/* subroutine for the calculation of the experimental unified glare rating (Fisekis 2003)*/
1244	+	float get_ugr_exp(pict * p, double lum_backg, double lum_a, int igs, int posindex_2)
1245	+	{
1246	+	float ugr_exp;
1247	+	double sum_glare;
1248	+	int i, i_glare;
1249
1250
1251	+	sum_glare = 0;
1252	+	i_glare = 0;
1253	+	for (i = 0; i <= igs; i++) {
1254	+
1255	+	if (pict_get_npix(p, i) > 0) {
1256	+	i_glare = i_glare + 1;
1257	+	sum_glare +=
1258	+	pow(1 /get_posindex(p, pict_get_av_posx(p, i), pict_get_av_posy(p, i), posindex_2),2) *pict_get_av_lum(p, i)* pict_get_av_omega(p, i);
1259	+
1260	+	}
1261	+	}
1262	+	ugr_exp =8 * log10(lum_a) + 8 * log10(sum_glare/lum_backg);
1263	+
1264	+	return ugr_exp;
1265	+
1266	+	}
1267	+	/* subroutine for the calculation of the unified glare probability (Hirning 2014)*/
1268	+	float get_ugp(pict * p, double lum_backg, int igs, int posindex_2)
1269	+	{
1270	+	float ugp;
1271	+	double sum_glare;
1272	+	int i, i_glare;
1273	+
1274	+
1275	+	sum_glare = 0;
1276	+	i_glare = 0;
1277	+	for (i = 0; i <= igs; i++) {
1278	+
1279	+	if (pict_get_npix(p, i) > 0) {
1280	+	i_glare = i_glare + 1;
1281	+	sum_glare +=
1282	+	pow(pict_get_av_lum(p, i) /
1283	+	get_posindex(p, pict_get_av_posx(p, i),
1284	+	pict_get_av_posy(p, i), posindex_2),
1285	+	2) * pict_get_av_omega(p, i);
1286	+
1287	+	}
1288	+	}
1289	+	ugp = 0.26 * log10(0.25 / lum_backg * sum_glare);
1290	+
1291	+	return ugp;
1292	+
1293	+	}
1294	+
1295	+
1296		/* subroutine for the calculation of the disability glare according to Poynter */
1297		float get_disability(pict * p, double lum_backg, int igs)
1298		{
#	Line 1218 | Line 1340 | float get_disability(pict * p, double lum_backg, int i
1340
1341		/* subroutine for the calculation of the cie glare index */
1342
1343	<	float
1222	<	get_cgi(pict * p, double E_v, double E_v_dir, int igs, int posindex_2)
1343	>	float get_cgi(pict * p, double E_v, double E_v_dir, int igs, int posindex_2)
1344		{
1345		float cgi;
1346		double sum_glare;
#	Line 1243 | Line 1364 | get_cgi(pict * p, double E_v, double E_v_dir, int igs,
1364		cgi = 8 * log10((2 * (1 + E_v_dir / 500) / E_v) * sum_glare);
1365
1366		return cgi;
1367	+	}
1368
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_con;
1375	+	double Lb;
1376	+
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	+
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	+
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! */
1389	+	float get_pgsv_sat(double E_v)
1390	+	{
1391	+	float pgsv_sat;
1392
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 1258 | Line 1433 | int main(int argc, char **argv)
1433		{
1434		#define CLINEMAX 4095 /* memory allocated for command line string */
1435		pict *p = pict_create();
1436	<	int     skip_second_scan,calcfast,age_corr,cut_view,cut_view_type,calc_vill, output, detail_out2, y1, fill, yfillmax, yfillmin,
1437	<	ext_vill, set_lum_max, set_lum_max2, task_color, i_splitstart,
1438	<	i_split, posindex_2, task_lum, checkfile, rval, i, i_max, x, y,x2,y2,
1439	<	igs, actual_igs, icol, xt, yt, change, before_igs, sgs, splithigh,uniform_gs,
1440	<	detail_out, posindex_picture, non_cos_lb, rx, ry, rmx, rmy,apply_disability,band_calc,band_color;
1441	<	double  lum_total_max,age_corr_factor,age,dgp_ext,dgp,low_light_corr,omega_cos_contr, setvalue, lum_ideal, E_v_contr, sigma,
1442	<	E_vl_ext, lum_max, new_lum_max, r_center,
1443	<	search_pix, a1, a2, a3, a4, a5, c3, c1, c2, r_split, max_angle,
1444	<	omegat, sang, E_v, E_v2, E_v_dir, avlum, act_lum, ang,
1445	<	l_max, lum_backg, lum_backg_cos, omega_sources, lum_sources,
1446	<	lum, lum_source,teta,Lveil_cie,Lveil_cie_sum,disability_thresh,u_r,u_g,u_b,band_angle,band_avlum;
1447	<	float lum_task, lum_thres, dgi,  vcp, cgi, ugr, limit,
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, 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_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,
1448	>	lum_mask[1],lum_mask_av,lum_mask_std[1],lum_mask_median[1],omega_mask,bbox[2],
1449	>	lum_band[1],lum_band_av,lum_band_std[1],lum_band_median[1],omega_band,bbox_band[2],
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,pgsv ,pgsv_sat,pgsv_con,
1454		abs_max, Lveil;
1455	<	char file_out[500], file_out2[500], version[500];
1455	>	char maskfile[500],file_out[500], file_out2[500], version[500];
1456		char *cline;
1457		VIEW userview = STDVIEW;
1458		int gotuserview = 0;
1459	+	struct muc_rvar* s_mask;
1460	+	s_mask = muc_rvar_create();
1461	+	muc_rvar_set_dim(s_mask, 1);
1462	+	muc_rvar_clear(s_mask);
1463	+	struct muc_rvar* s_band;
1464	+	s_band = muc_rvar_create();
1465	+	muc_rvar_set_dim(s_band, 1);
1466	+	muc_rvar_clear(s_band);
1467	+	struct muc_rvar* s_z1;
1468	+	s_z1 = muc_rvar_create();
1469	+	muc_rvar_set_dim(s_z1, 1);
1470	+	muc_rvar_clear(s_z1);
1471
1472	+	struct muc_rvar* s_z2;
1473	+	s_z2 = muc_rvar_create();
1474	+	muc_rvar_set_dim(s_z2, 1);
1475	+	muc_rvar_clear(s_z2);
1476	+
1477	+	struct muc_rvar* s_noposweight;
1478	+	s_noposweight = muc_rvar_create();
1479	+	muc_rvar_set_dim(s_noposweight, 1);
1480	+	muc_rvar_clear(s_noposweight);
1481	+
1482	+	struct muc_rvar* s_posweight;
1483	+	s_posweight = muc_rvar_create();
1484	+	muc_rvar_set_dim(s_posweight, 1);
1485	+	muc_rvar_clear(s_posweight);
1486	+
1487	+	struct muc_rvar* s_posweight2;
1488	+	s_posweight2 = muc_rvar_create();
1489	+	muc_rvar_set_dim(s_posweight2, 1);
1490	+	muc_rvar_clear(s_posweight2);
1491	+
1492		/*set required user view parameters to invalid values*/
1493	<	uniform_gs = 0;
1493	>	dir_ill=0.0;
1494	>	delta_E=0.0;
1495	>	no_glaresources=0;
1496	>	n_corner_px=0;
1497	>	zero_corner_px=0;
1498	>	force=0;
1499	>	dist=0.0;
1500	>	u_r=0.33333;
1501	>	u_g=0.33333;
1502	>	u_b=0.33333;
1503	>	band_angle=0;
1504	>	angle_z1=0;
1505	>	angle_z2=0;
1506	>	x_zone=0;
1507	>	y_zone=0;
1508	>	per_75_z2=0;
1509	>	per_95_z2=0;
1510	>	lum_pos_mean=0;
1511	>	lum_pos2_mean=0;
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;
1519	>	omega_z1 = 0.0;
1520	>	lum_z2_av = 0.0;
1521	>	omega_z2 = 0.0 ;
1522	>	i_z1 = 0;
1523	>	i_z2 = 0;
1524	>	zones = 0;
1525	>	lum_z2_av = 0.0;
1526	>	uniform_gs = 0;
1527		apply_disability = 0;
1528		disability_thresh = 0;
1529		Lveil_cie_sum=0.0;
#	Line 1297 | Line 1543 | int main(int argc, char **argv)
1543		omegat = 0.0;
1544		yt = 0;
1545		xt = 0;
1546	+	x_disk=0;
1547	+	y_disk=0;
1548	+	angle_disk=0;
1549		yfillmin = 0;
1550		yfillmax = 0;
1551		cut_view = 0;
#	Line 1305 | 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 1323 | 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;
1580		set_lum_max = 0;
1581		set_lum_max2 = 0;
1582	+	img_corr=0;
1583		abs_max = 0;
1584		progname = argv[0];
1585		E_v_contr = 0.0;
1586	<	strcpy(version, "1.15 release 05.02.2015 by J.Wienold");
1586	>	strcpy(version, "1.19 release 09.12.2015 by J.Wienold");
1587		low_light_corr=1.0;
1588		output = 0;
1589		calc_vill = 0;
1590		band_avlum = -99;
1591		band_calc = 0;
1592	+	masking = 0;
1593	+	lum_mask[0]=0.0;
1594	+	lum_mask_av=0.0;
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 1379 | Line 1638 | int main(int argc, char **argv)
1638		printf("age factor not supported any more \n");
1639		exit(1);
1640		break;
1641	+	case 'A':
1642	+	masking = 1;
1643	+	detail_out = 1;
1644	+	strcpy(maskfile, argv[++i]);
1645	+	pict_read(pm,maskfile );
1646	+
1647	+	break;
1648		case 'b':
1649		lum_thres = atof(argv[++i]);
1650	+	thres_activate = 1;
1651		break;
1652		case 'c':
1653		checkfile = 1;
#	Line 1398 | Line 1665 | int main(int argc, char **argv)
1665		case 's':
1666		sgs = 1;
1667		break;
1668	+	case 'f':
1669	+	force = 1;
1670	+	break;
1671		case 'k':
1672		apply_disability = 1;
1673		disability_thresh = atof(argv[++i]);
#	Line 1406 | Line 1676 | int main(int argc, char **argv)
1676		case 'p':
1677		posindex_picture = 1;
1678		break;
1679	+	case 'P':
1680	+	posindex_2 = atoi(argv[++i]);
1681	+	break;
1682	+
1683
1684		case 'y':
1685		splithigh = 1;
#	Line 1436 | Line 1710 | int main(int argc, char **argv)
1710		detail_out2 = 1;
1711		break;
1712		case 'm':
1713	+	img_corr = 1;
1714		set_lum_max = 1;
1715		lum_max = atof(argv[++i]);
1716		new_lum_max = atof(argv[++i]);
#	Line 1444 | Line 1719 | int main(int argc, char **argv)
1719		break;
1720
1721		case 'M':
1722	+	img_corr = 1;
1723		set_lum_max2 = 1;
1724		lum_max = atof(argv[++i]);
1725		E_vl_ext = atof(argv[++i]);
1726		strcpy(file_out2, argv[++i]);
1727		/*                      printf("max lum set to %f\n",new_lum_max);*/
1728		break;
1729	<	case 'n':
1729	>	case 'N':
1730	>	img_corr = 1;
1731	>	set_lum_max2 = 2;
1732	>	x_disk = atoi(argv[++i]);
1733	>	y_disk = atoi(argv[++i]);
1734	>	angle_disk = atof(argv[++i]);
1735	>	E_vl_ext = atof(argv[++i]);
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	>
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 1469 | Line 1771 | int main(int argc, char **argv)
1771		omegat = atof(argv[++i]);
1772		task_color = 1;
1773		break;
1774	+	case 'l':
1775	+	zones = 1;
1776	+	x_zone = atoi(argv[++i]);
1777	+	y_zone = atoi(argv[++i]);
1778	+	angle_z1 = atof(argv[++i]);
1779	+	angle_z2 = -1;
1780	+	break;
1781	+	case 'L':
1782	+	zones = 2;
1783	+	x_zone = atoi(argv[++i]);
1784	+	y_zone = atoi(argv[++i]);
1785	+	angle_z1 = atof(argv[++i]);
1786	+	angle_z2 = atof(argv[++i]);
1787	+	break;
1788	+
1789	+
1790		case 'B':
1791		band_calc = 1;
1792		band_angle = atof(argv[++i]);
#	Line 1505 | 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 1534 | 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) {
1877	+	fprintf(stderr,  " masking and zoning cannot be activated at the same time!\n");
1878	+	exit(1);
1879	+	}
1880	+
1881		/* read picture file */
1882		if (i == argc) {
1883		SET_FILE_BINARY(stdin);
#	Line 1573 | Line 1911 | if (output == 1 && ext_vill == 1) {
1911		return EXIT_FAILURE;
1912		}
1913
1914	+
1915	+
1916		/*                fprintf(stderr,"Picture read!");*/
1917
1918		/* command line for output picture*/
1919
1920		pict_set_comment(p, cline);
1921
1922	+	/* several checks */
1923		if (p->view.type == VT_PAR) {
1924	<	fprintf(stderr, "wrong view type! must not be parallel ! \n");
1924	>	fprintf(stderr, "error: wrong view type! must not be parallel ! \n");
1925		exit(1);
1926		}
1927
1928	+	if (p->view.type == VT_PER) {
1929	+	fprintf(stderr, "warning: image has perspective view type specified in header ! \n");
1930	+	}
1931
1932	+	if (masking == 1) {
1933	+
1934	+	if (!pict_get_xsize(p)==pict_get_xsize(pm) || !pict_get_ysize(p)==pict_get_ysize(pm)) {
1935	+	fprintf(stderr, "error: masking image has other resolution than main image ! \n");
1936	+	fprintf(stderr, "size must be identical \n");
1937	+	printf("resolution main image : %dx%d\n",pict_get_xsize(p),pict_get_ysize(p));
1938	+	printf("resolution masking image : %dx%d\n",pict_get_xsize(pm),pict_get_ysize(pm));
1939	+	exit(1);
1940	+
1941	+	}
1942	+
1943	+	}
1944	+	/*      printf("resolution: %dx%d\n",pict_get_xsize(p),pict_get_ysize(p)); */
1945	+
1946		/* Check size of search radius */
1947		rmx = (pict_get_xsize(p) / 2);
1948		rmy = (pict_get_ysize(p) / 2);
#	Line 1629 | Line 1987 | if (output == 1 && ext_vill == 1) {
1987		avlum = 0.0;
1988		pict_new_gli(p);
1989		igs = 0;
1990	+	pict_get_z1_gsn(p,igs) = 0;
1991	+	pict_get_z2_gsn(p,igs) = 0;
1992
1993	+
1994		/* cut out GUTH field of view and exit without glare evaluation */
1995		if (cut_view==2) {
1996		if (cut_view_type==1) {
#	Line 1696 | Line 2057 | if (cut_view==2) {
2057		if (calcfast == 0) {
2058		for (x = 0; x < pict_get_xsize(p); x++)
2059		for (y = 0; y < pict_get_ysize(p); y++) {
2060	+	lum = luminance(pict_get_color(p, x, y));
2061	+	dist=sqrt((x-rmx)*(x-rmx)+(y-rmy)*(y-rmy));
2062	+	if (dist>((rmx+rmy)/2)) {
2063	+	n_corner_px=n_corner_px+1;
2064	+	if (lum<7.0) {
2065	+	zero_corner_px=zero_corner_px+1;
2066	+	}
2067	+	}
2068	+
2069	+
2070		if (pict_get_hangle(p, x, y, p->view.vdir, p->view.vup, &ang)) {
2071		if (pict_is_validpixel(p, x, y)) {
1701	–	lum = luminance(pict_get_color(p, x, y));
2072		if (fill == 1 && y >= yfillmax) {
2073		y1 = y - 1;
2074		lum = luminance(pict_get_color(p, x, y1));
#	Line 1711 | Line 2081 | if (cut_view==2) {
2081		abs_max = lum;
2082		}
2083		/* set luminance restriction, if -m is set */
2084	<	if (set_lum_max == 1 && lum >= lum_max) {
2085	<	pict_get_color(p, x, y)[RED] = new_lum_max / 179.0;
2086	<	pict_get_color(p, x, y)[GRN] = new_lum_max / 179.0;
2087	<	pict_get_color(p, x, y)[BLU] = new_lum_max / 179.0;
2088	<	/*                                    printf("old lum, new lum %f %f \n",lum,luminance(pict_get_color(p,x,y))); */
2089	<	lum = luminance(pict_get_color(p, x, y));
2084	>	if (img_corr == 1 ) {
2085	>	if (set_lum_max == 1 && lum >= lum_max) {
2086	>	pict_get_color(p, x, y)[RED] = new_lum_max / 179.0;
2087	>	pict_get_color(p, x, y)[GRN] = new_lum_max / 179.0;
2088	>	pict_get_color(p, x, y)[BLU] = new_lum_max / 179.0;
2089	>	/*                                                      printf("old lum, new lum %f %f \n",lum,luminance(pict_get_color(p,x,y))); */
2090	>	lum = luminance(pict_get_color(p, x, y));
2091	>	}
2092	>	if (set_lum_max2 == 1 && lum >= lum_max) {
2093
2094	<	}
2095	<	if (set_lum_max2 == 1 && lum >= lum_max) {
2094	>	E_v_contr +=DOT(p->view.vdir,pict_get_cached_dir(p, x,y)) * pict_get_omega(p,x,y)* lum;
2095	>	omega_cos_contr +=DOT(p->view.vdir,pict_get_cached_dir(p, x,y)) * pict_get_omega(p,x,y)* 1;
2096	>	}
2097	>	if (set_lum_max2 == 2 ) {
2098	>	r_actual =acos(DOT(pict_get_cached_dir(p, x_disk, y_disk), pict_get_cached_dir(p, x, y))) * 2;
2099	>	if (x_disk == x && y_disk==y ) r_actual=0.0;
2100
2101	<	E_v_contr +=
2102	<	DOT(p->view.vdir,
2103	<	pict_get_cached_dir(p, x,
2104	<	y)) * pict_get_omega(p,
2105	<	x,
2106	<	y)
2107	<	* lum;
2108	<	omega_cos_contr +=
2109	<	DOT(p->view.vdir,
2110	<	pict_get_cached_dir(p, x,
1734	<	y)) * pict_get_omega(p,
1735	<	x,
1736	<	y)
1737	<	* 1;
2101	>	act_lum = luminance(pict_get_color(p, x, y));
2102	>
2103	>	if (r_actual <= angle_disk) {
2104	>	E_v_contr +=DOT(p->view.vdir,pict_get_cached_dir(p, x,y)) * pict_get_omega(p,x,y)* lum;
2105	>	omega_cos_contr += DOT(p->view.vdir,pict_get_cached_dir(p, x,y)) * pict_get_omega(p,x,y)* 1;
2106	>	}
2107	>
2108	>
2109	>
2110	>	}
2111		}
2112	+	om = pict_get_omega(p, x, y);
2113	+	sang += om;
2114	+	E_v += DOT(p->view.vdir, pict_get_cached_dir(p, x,y)) * om *lum;
2115	+	avlum += om * lum;
2116	+	pos=get_posindex(p, x, y,posindex_2);
2117
2118	+	lum_pos[0]=lum;
2119	+	muc_rvar_add_sample(s_noposweight,lum_pos);
2120	+	lum_pos[0]=lum/pos;
2121	+	lum_pos_mean+=lum_pos[0]*om;
2122	+	muc_rvar_add_sample(s_posweight, lum_pos);
2123	+	lum_pos[0]=lum_pos[0]/pos;
2124	+	lum_pos2_mean+=lum_pos[0]*om;
2125	+	muc_rvar_add_sample(s_posweight2,lum_pos);
2126
2127	<	sang += pict_get_omega(p, x, y);
2128	<	E_v +=
1743	<	DOT(p->view.vdir,
1744	<	pict_get_cached_dir(p, x,
1745	<	y)) * pict_get_omega(p, x,
1746	<	y) *
1747	<	lum;
1748	<	avlum +=
1749	<	pict_get_omega(p, x,
1750	<	y) * luminance(pict_get_color(p, x,
1751	<	y));
2127	>
2128	>
2129		} else {
2130		pict_get_color(p, x, y)[RED] = 0.0;
2131		pict_get_color(p, x, y)[GRN] = 0.0;
2132		pict_get_color(p, x, y)[BLU] = 0.0;
2133
2134		}
2135	<	}
2135	>	}else {
2136	>	pict_get_color(p, x, y)[RED] = 0.0;
2137	>	pict_get_color(p, x, y)[GRN] = 0.0;
2138	>	pict_get_color(p, x, y)[BLU] = 0.0;
2139	>
2140	>	}
2141		}
2142
2143	<	if (set_lum_max2 == 1 && E_v_contr > 0 && (E_vl_ext - E_v) > 0) {
2143	>	/* check if image has black corners AND a perspective view */
2144
2145	+	if (zero_corner_px/n_corner_px > 0.70 && (p->view.type == VT_PER) ) {
2146	+	printf (" corner pixels are to  %f %% black! \n",zero_corner_px/n_corner_px*100);
2147	+	printf("error! The image has black corners AND header contains a perspective view type definition !!!\n") ;
2148	+
2149	+	if (force==0) {
2150	+	printf("stopping...!!!!\n") ;
2151	+
2152	+	exit(1);
2153	+	}
2154	+	}
2155	+
2156	+
2157	+
2158	+
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 ) || 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") ;
2168	+	}
2169	+
2170		if (lum_ideal > 0 && lum_ideal < setvalue) {
2171		setvalue = lum_ideal;
2172		}
2173	<	printf
1768	<	("change luminance values!! lum_ideal,setvalue,E_vl_ext,E_v,E_v_contr %f  %f %f %f %f\n",
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	<
2178	>
2179		for (x = 0; x < pict_get_xsize(p); x++)
2180		for (y = 0; y < pict_get_ysize(p); y++) {
2181		if (pict_get_hangle
2182		(p, x, y, p->view.vdir, p->view.vup, &ang)) {
2183		if (pict_is_validpixel(p, x, y)) {
2184		lum = luminance(pict_get_color(p, x, y));
2185	+
2186	+
2187		if (set_lum_max2 == 1 && lum >= lum_max) {
2188
2189		pict_get_color(p, x, y)[RED] =
#	Line 1784 | 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 >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	+
2200	+	if (r_actual <= angle_disk) {
2201	+
2202	+	pict_get_color(p, x, y)[RED] =
2203	+	setvalue / 179.0;
2204	+	pict_get_color(p, x, y)[GRN] =
2205	+	setvalue / 179.0;
2206	+	pict_get_color(p, x, y)[BLU] =
2207	+	setvalue / 179.0;
2208	+
2209	+	}
2210		}
2211	+	}
2212		}
2213		}
2214		}
2215	+
2216
2217		pict_write(p, file_out2);
2218	<
2218	>	exit(1);
2219		}
2220		if (set_lum_max == 1) {
2221		pict_write(p, file_out2);
#	Line 1850 | Line 2275 | if (cut_view==1) {
2275		lum_source = lum_thres;
2276		}
2277
2278	<	/*     printf("Ev, avlum, Omega_tot%f %f %f \n",E_v, avlum, sang); */
2278	>	/*     printf("Ev, avlum, Omega_tot%f %f %f \n",E_v, avlum, sang); */
2279
2280		/* first glare source scan: find primary glare sources */
2281	<	for (x = 0; x < pict_get_xsize(p); x++)
2281	>	for (x = 0; x < pict_get_xsize(p); x++) {
2282	>	/*                lastpixelwas_gs=0; */
2283		for (y = 0; y < pict_get_ysize(p); y++) {
2284		if (pict_get_hangle(p, x, y, p->view.vdir, p->view.vup, &ang)) {
2285		if (pict_is_validpixel(p, x, y)) {
#	Line 1862 | Line 2288 | if (cut_view==1) {
2288		if (act_lum >lum_total_max) {
2289		lum_total_max=act_lum;
2290		}
2291	<
2292	<	actual_igs =
2293	<	find_near_pgs(p, x, y, max_angle, 0, igs, 1);
2291	>	/* speed improvement first scan: when last pixel was glare source, then use this glare source number instead of searching.
2292	>	has been deactivated with v1.25
2293	>
2294	>	if (lastpixelwas_gs==0 || search_pix <= 1.0 ) { */
2295	>	actual_igs = find_near_pgs(p, x, y, max_angle, 0, igs, 1);
2296	>	/*  }*/
2297		if (actual_igs == 0) {
2298		igs = igs + 1;
2299		pict_new_gli(p);
2300		pict_get_lum_min(p, igs) = HUGE_VAL;
2301		pict_get_Eglare(p,igs) = 0.0;
2302		pict_get_Dglare(p,igs) = 0.0;
2303	+	pict_get_z1_gsn(p,igs) = 0;
2304	+	pict_get_z2_gsn(p,igs) = 0;
2305		actual_igs = igs;
2306	+
2307		}
2308	<	icol = setglcolor(p, x, y, actual_igs, uniform_gs, u_r, u_g , u_b);
2308	>	/* no coloring any more here            icol = setglcolor(p, x, y, actual_igs, uniform_gs, u_r, u_g , u_b); */
2309		pict_get_gsn(p, x, y) = actual_igs;
2310		pict_get_pgs(p, x, y) = 1;
2311		add_pixel_to_gs(p, x, y, actual_igs);
2312	+	/*                                                lastpixelwas_gs=actual_igs; */
2313
2314		} else {
2315		pict_get_pgs(p, x, y) = 0;
2316	+	/*                                               lastpixelwas_gs=0;*/
2317		}
2318		}
2319		}
2320	<	}
2320	>	}
2321	>	}
2322		/*      pict_write(p,"firstscan.pic");   */
2323
2324	<	if (calcfast == 1 ) {
2324	>
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;
2334	+	i = 0;
2335		while (change == 1 && skip_second_scan==0) {
2336		change = 0;
2337	<	for (x = 0; x < pict_get_xsize(p); x++)
2337	>	i = i+1;
2338	>	for (x = 0; x < pict_get_xsize(p); x++) {
2339		for (y = 0; y < pict_get_ysize(p); y++) {
1899	–	before_igs = pict_get_gsn(p, x, y);
2340		if (pict_get_hangle
2341		(p, x, y, p->view.vdir, p->view.vup, &ang)) {
2342	<	if (pict_is_validpixel(p, x, y) && before_igs > 0) {
2342	>	checkpixels=1;
2343	>	before_igs = pict_get_gsn(p, x, y);
2344	>
2345	>	/* speed improvement: search for other nearby glare sources only, when at least one adjacend pixel has another glare source number  */
2346	>	if (x >1 && x<pict_get_xsize(p)-2 && y >1 && y<pict_get_ysize(p)-2 )   {
2347	>	x1=x-1;
2348	>	x2=x+1;
2349	>	y1=y-1;
2350	>	y2=y+1;
2351	>	if (before_igs == pict_get_gsn(p, x1, y) && before_igs == pict_get_gsn(p, x2, y) && before_igs == pict_get_gsn(p, x, y1) && before_igs == pict_get_gsn(p, x, y2) && before_igs == pict_get_gsn(p, x1, y1) && before_igs == pict_get_gsn(p, x2, y1) && before_igs == pict_get_gsn(p, x1, y2) && before_igs == pict_get_gsn(p, x2, y2) ) {
2352	>	checkpixels = 0;
2353	>	actual_igs = before_igs;
2354	>	}  }
2355	>
2356	>
2357	>	if (pict_is_validpixel(p, x, y) && before_igs > 0 && checkpixels==1) {
2358		actual_igs =
2359		find_near_pgs(p, x, y, max_angle, before_igs,
2360		igs, 1);
#	Line 1908 | Line 2363 | if (calcfast == 1 ) {
2363		}
2364		if (before_igs > 0) {
2365		actual_igs = pict_get_gsn(p, x, y);
2366	<	icol = setglcolor(p, x, y, actual_igs, uniform_gs, u_r, u_g , u_b);
2366	>	/*      icol = setglcolor(p, x, y, actual_igs, uniform_gs, u_r, u_g , u_b); */
2367		}
2368
2369		}
2370		}
2371		}
2372		/*      pict_write(p,"secondscan.pic");*/
2373	+	}}
2374
1919	–	}
1920	–
2375		/*smoothing: add secondary glare sources */
2376		i_max = igs;
2377		if (sgs == 1) {
#	Line 1974 | Line 2428 | if (calcfast == 1 ) {
2428		if (i_splitstart == (igs + 1)) {
2429		igs = igs + 1;
2430		pict_new_gli(p);
2431	+	pict_get_z1_gsn(p,igs) = 0;
2432	+	pict_get_z2_gsn(p,igs) = 0;
2433	+
2434		pict_get_Eglare(p,igs) = 0.0;
2435		pict_get_Dglare(p,igs) = 0.0;
2436		pict_get_lum_min(p, igs) =
#	Line 1987 | Line 2444 | if (calcfast == 1 ) {
2444		if (i_split == 0) {
2445		igs = igs + 1;
2446		pict_new_gli(p);
2447	+	pict_get_z1_gsn(p,igs) = 0;
2448	+	pict_get_z2_gsn(p,igs) = 0;
2449	+
2450		pict_get_Eglare(p,igs) = 0.0;
2451		pict_get_Dglare(p,igs) = 0.0;
2452		pict_get_lum_min(p, igs) =
#	Line 2023 | Line 2483 | if (calcfast == 1 ) {
2483		if (before_igs > 0) {
2484		actual_igs =
2485		pict_get_gsn(p, x, y);
2486	<	icol =
2486	>	/*                                                                                     icol =
2487		setglcolor(p, x, y,
2488	<	actual_igs, uniform_gs, u_r, u_g , u_b);
2488	>	actual_igs, uniform_gs, u_r, u_g , u_b);*/
2489		}
2490
2491		}
#	Line 2040 | Line 2500 | if (calcfast == 1 ) {
2500		}
2501		}
2502
2503	<	/* calculation of direct vertical illuminance for CGI and for disability glare*/
2503	>	/* calculation of direct vertical illuminance for CGI and for disability glare, coloring glare sources*/
2504
2505	<
2046	<	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)) {
2509		if (pict_is_validpixel(p, x, y)) {
2510		if (pict_get_gsn(p, x, y) > 0) {
2511	<	pict_get_Eglare(p,pict_get_gsn(p, x, y)) +=
2512	<	DOT(p->view.vdir, pict_get_cached_dir(p, x, y))
2513	<	* pict_get_omega(p, x, y)
2514	<	* luminance(pict_get_color(p, x, y));
2515	<	E_v_dir +=
2057	<	DOT(p->view.vdir, pict_get_cached_dir(p, x, y))
2058	<	* pict_get_omega(p, x, y)
2059	<	* luminance(pict_get_color(p, x, y));
2511	>	actual_igs = pict_get_gsn(p, x, y);
2512	>	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));
2513	>	pict_get_Eglare(p,actual_igs ) = pict_get_Eglare(p,actual_igs ) + delta_E;
2514	>	E_v_dir = E_v_dir + delta_E;
2515	>	setglcolor(p, x, y, actual_igs, uniform_gs, u_r, u_g , u_b);
2516		}
2517		}
2518		}
2519		}
2520		lum_backg_cos = (E_v - E_v_dir) / 3.1415927;
2521	<	lum_backg = lum_backg_cos;
2521	>
2522		}
2523	<	/* calc of band luminance if applied */
2523	>	/* calc of band luminance distribution if applied */
2524		if (band_calc == 1) {
2525	<	band_avlum=get_band_lum( p, band_angle, band_color);
2526	<	}
2525	>	x_max = pict_get_xsize(p) - 1;
2526	>	y_max = pict_get_ysize(p) - 1;
2527	>	y_mid = (int)(y_max/2);
2528	>	for (x = 0; x <= x_max; x++)
2529	>	for (y = 0; y <= y_max; y++) {
2530	>	if (pict_get_hangle(p, x, y, p->view.vdir, p->view.vup, &ang)) {
2531	>	if (pict_is_validpixel(p, x, y)) {
2532
2533	+	r_actual = acos(DOT(pict_get_cached_dir(p, x, y_mid), pict_get_cached_dir(p, x, y))) ;
2534	+	act_lum = luminance(pict_get_color(p, x, y));
2535	+
2536	+	if ((r_actual <= band_angle) || (y == y_mid) ) {
2537	+	lum_band[0]=luminance(pict_get_color(p, x, y));
2538	+	muc_rvar_add_sample(s_band, lum_band);
2539	+	act_lum = luminance(pict_get_color(p, x, y));
2540	+	lum_band_av += pict_get_omega(p, x, y) * act_lum;
2541	+	omega_band += pict_get_omega(p, x, y);
2542	+	if (band_color == 1) {
2543	+	pict_get_color(p, x, y)[RED] = 0.0;
2544	+	pict_get_color(p, x, y)[GRN] = act_lum / WHTEFFICACY / CIE_gf;
2545	+	pict_get_color(p, x, y)[BLU] = 0.0;
2546	+	}
2547	+	}
2548	+	}}}
2549	+	lum_band_av=lum_band_av/omega_band;
2550	+	muc_rvar_get_vx(s_band,lum_band_std);
2551	+	muc_rvar_get_percentile(s_band,lum_band_median,0.75);
2552	+	per_75_band=lum_band_median[0];
2553	+	muc_rvar_get_percentile(s_band,lum_band_median,0.95);
2554	+	per_95_band=lum_band_median[0];
2555	+	muc_rvar_get_median(s_band,lum_band_median);
2556	+	muc_rvar_get_bounding_box(s_band,bbox_band);
2557	+
2558	+	printf ("band:band_omega,band_av_lum,band_median_lum,band_std_lum,band_perc_75,band_perc_95,band_lum_min,band_lum_max: %f %f %f %f %f %f %f %f\n",omega_band,lum_band_av,lum_band_median[0],sqrt(lum_band_std[0]),per_75_band,per_95_band,bbox_band[0],bbox_band[1] );
2559	+	/*      av_lum = av_lum / omega_sum; */
2560	+	/*    printf("average luminance of band %f \n",av_lum);*/
2561	+	/*      printf("solid angle of band %f \n",omega_sum);*/
2562	+	}
2563	+
2564		/*printf("total number of glare sources: %i\n",igs); */
2565		lum_sources = 0;
2566		omega_sources = 0;
2567	+	i=0;
2568		for (x = 0; x <= igs; x++) {
2569	+	if (pict_get_npix(p, x) > 0) {
2570		lum_sources += pict_get_av_lum(p, x) * pict_get_av_omega(p, x);
2571		omega_sources += pict_get_av_omega(p, x);
2572	+	i=i+1;
2573	+	}}
2574	+
2575	+	if (sang == omega_sources) {
2576	+	lum_backg =avlum;
2577	+	} else {
2578	+	lum_backg =(sang * avlum - lum_sources) / (sang - omega_sources);
2579		}
2580	<	if (non_cos_lb == 1) {
2581	<	lum_backg =
2582	<	(sang * avlum - lum_sources) / (sang - omega_sources);
2580	>
2581	>
2582	>	if (i == 0) {
2583	>	lum_sources =0.0;
2584	>	} else { lum_sources=lum_sources/omega_sources;
2585		}
2586	+
2587	+
2588	+
2589	+	if (non_cos_lb == 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);
2601	+	}
2602	+	*/
2603	+
2604		/* print detailed output */
2605	+
2606	+
2607	+
2608	+	/* masking */
2609	+
2610	+	if ( masking == 1 ) {
2611	+
2612	+
2613	+	for (x = 0; x < pict_get_xsize(p); x++)
2614	+	for (y = 0; y < pict_get_ysize(p); y++) {
2615	+	if (pict_get_hangle(p, x, y, p->view.vdir, p->view.vup, &ang)) {
2616	+	if (pict_is_validpixel(p, x, y)) {
2617	+	if (luminance(pict_get_color(pm, x, y))>0.1) {
2618	+	/*                                         printf ("hallo %i %i %f \n",x,y,luminance(pict_get_color(pm, x, y)));*/
2619	+	i_mask=i_mask+1;
2620	+	lum_mask[0]=luminance(pict_get_color(p, x, y));
2621	+	/* printf ("%f \n",lum_mask[0]);*/
2622	+	muc_rvar_add_sample(s_mask, lum_mask);
2623	+	omega_mask += pict_get_omega(p, x, y);
2624	+	lum_mask_av += pict_get_omega(p, x, y)* luminance(pict_get_color(p, x, y));
2625	+	E_v_mask +=DOT(p->view.vdir, pict_get_cached_dir(p, x, y))*pict_get_omega(p, x, y)* luminance(pict_get_color(p, x, y));
2626	+
2627	+	pict_get_color(pm, x, y)[RED] = luminance(pict_get_color(p, x, y))/179.0;
2628	+	pict_get_color(pm, x, y)[GRN] = luminance(pict_get_color(p, x, y))/179.0;
2629	+	pict_get_color(pm, x, y)[BLU] = luminance(pict_get_color(p, x, y))/179.0;
2630	+
2631	+	} else {
2632	+	pict_get_color(p, x, y)[RED] = 0.0 ;
2633	+	pict_get_color(p, x, y)[GRN] = 0.0 ;
2634	+	pict_get_color(p, x, y)[BLU] = 0.0 ;
2635	+	}
2636	+	}
2637	+	}
2638	+	}
2639	+	/* Average luminance in masking area (weighted by solid angle) */
2640	+	lum_mask_av=lum_mask_av/omega_mask;
2641	+	muc_rvar_get_vx(s_mask,lum_mask_std);
2642	+	muc_rvar_get_percentile(s_mask,lum_mask_median,0.75);
2643	+	per_75_mask=lum_mask_median[0];
2644	+	muc_rvar_get_percentile(s_mask,lum_mask_median,0.95);
2645	+	per_95_mask=lum_mask_median[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	+
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);
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	+
2667	+	}
2668	+
2669	+	/* zones */
2670	+
2671	+	if ( zones > 0 ) {
2672	+
2673	+	for (x = 0; x < pict_get_xsize(p); x++)
2674	+	for (y = 0; y < pict_get_ysize(p); y++) {
2675	+	if (pict_get_hangle(p, x, y, p->view.vdir, p->view.vup, &ang)) {
2676	+	if (pict_is_validpixel(p, x, y)) {
2677	+
2678	+
2679	+	r_actual =acos(DOT(pict_get_cached_dir(p, x, y), pict_get_cached_dir(p, x_zone, y_zone))) * 2;
2680	+	lum_actual = luminance(pict_get_color(p, x, y));
2681	+	act_gsn=pict_get_gsn(p, x, y);
2682	+	/*     printf("1act_gsn,pict_get_z1_gsn,pict_get_z2_gsn %i %f %f \n",act_gsn,pict_get_z1_gsn(p,act_gsn),pict_get_z2_gsn(p,act_gsn));*/
2683	+
2684	+	/*zone1 */
2685	+	if (r_actual <= angle_z1) {
2686	+	i_z1=i_z1+1;
2687	+	lum_z1[0]=lum_actual;
2688	+	muc_rvar_add_sample(s_z1, lum_z1);
2689	+	omega_z1 += pict_get_omega(p, x, y);
2690	+	lum_z1_av += pict_get_omega(p, x, y)* lum_actual;
2691	+	setglcolor(p,x,y,1,1 , 0.66, 0.01 ,0.33);
2692	+	/*check if separation gsn already exist */
2693	+
2694	+	if (act_gsn > 0){
2695	+	if (pict_get_z1_gsn(p,act_gsn) == 0) {
2696	+	pict_new_gli(p);
2697	+	igs = igs + 1;
2698	+	pict_get_z1_gsn(p,act_gsn) = igs*1.000;
2699	+	pict_get_z1_gsn(p,igs) = -1.0;
2700	+	pict_get_z2_gsn(p,igs) = -1.0;
2701	+	splitgs=(int)(pict_get_z1_gsn(p,act_gsn));
2702	+	/*                                                                                        printf ("Glare source %i gets glare source %i in zone 1 : %i %f %f \n",act_gsn,igs,splitgs,pict_get_z1_gsn(p,act_gsn),pict_get_z1_gsn(p,igs)); */
2703	+	}
2704	+	splitgs=(int)(pict_get_z1_gsn(p,act_gsn));
2705	+	/*          printf("splitgs%i \n",splitgs);       */
2706	+	split_pixel_from_gs(p, x, y, splitgs, uniform_gs, u_r, u_g , u_b);
2707	+	}
2708	+	}
2709	+	/*zone2 */
2710	+
2711	+	if (r_actual > angle_z1 && r_actual<= angle_z2 ) {
2712	+
2713	+	i_z2=i_z2+1;
2714	+	lum_z2[0]=lum_actual;
2715	+	muc_rvar_add_sample(s_z2, lum_z2);
2716	+	omega_z2 +=   pict_get_omega(p, x, y);
2717	+	lum_z2_av += pict_get_omega(p, x, y)* lum_actual;
2718	+	setglcolor(p,x,y,1,1 , 0.65, 0.33 ,0.02);
2719	+	/*                                                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));
2720	+	*/                                                if (act_gsn > 0){
2721	+	if (pict_get_z2_gsn(p,act_gsn) == 0) {
2722	+	pict_new_gli(p);
2723	+	igs = igs + 1;
2724	+	pict_get_z2_gsn(p,act_gsn) = igs*1.000;
2725	+	pict_get_z1_gsn(p,igs) = -2.0;
2726	+	pict_get_z2_gsn(p,igs) = -2.0;
2727	+	/*                                                                                         printf ("Glare source %i gets glare source %i in zone 2 \n",act_gsn,igs ); */
2728	+	}
2729	+	splitgs=(int)(pict_get_z2_gsn(p,act_gsn));
2730	+	/*                                              printf("splitgs %i \n",splitgs);*/
2731	+	split_pixel_from_gs(p, x, y, splitgs, uniform_gs, u_r, u_g , u_b);
2732	+	}
2733	+	}
2734	+
2735	+	}
2736	+	} }
2737	+	/* Average luminance in zones (weighted by solid angle), calculate percentiles, median min and max in zones  */
2738	+	if (zones == 2 ) {
2739	+	lum_z2_av=lum_z2_av/omega_z2;
2740	+	muc_rvar_get_vx(s_z2,lum_z2_std);
2741	+	muc_rvar_get_percentile(s_z2,lum_z2_median,0.75);
2742	+	per_75_z2=lum_z2_median[0];
2743	+	muc_rvar_get_percentile(s_z2,lum_z2_median,0.95);
2744	+	per_95_z2=lum_z2_median[0];
2745	+	muc_rvar_get_median(s_z2,lum_z2_median);
2746	+	muc_rvar_get_bounding_box(s_z2,bbox_z2);
2747	+	}
2748	+	lum_z1_av=lum_z1_av/omega_z1;
2749	+	muc_rvar_get_vx(s_z1,lum_z1_std);
2750	+	muc_rvar_get_percentile(s_z1,lum_z1_median,0.75);
2751	+	per_75_z1=lum_z1_median[0];
2752	+	muc_rvar_get_percentile(s_z1,lum_z1_median,0.95);
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	+	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	+	} }
2764	+
2765	+	}
2766	+
2767		i = 0;
2768		for (x = 0; x <= igs; x++) {
2769		/* resorting glare source numbers */
#	Line 2089 | Line 2771 | if (calcfast == 1 ) {
2771		i = i + 1;
2772		}
2773		}
2774	+	no_glaresources=i;
2775
2776	+	/* glare sources */
2777	+	if (detail_out == 1) {
2778
2094	–
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 \n",
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);
2782		if (i == 0) {
2783	<	printf("%i %f %f %f %f %.10f %f %f %f %f %f %f\n", i, 0.0, 0.0,
2784	<	0.0, 0.0, 0.0, 0.0, lum_backg, lum_task, E_v, E_v_dir,
2101	<	abs_max);
2783	>	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,
2784	>	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);
2785
2786		} else {
2787		i = 0;
#	Line 2118 | Line 2801 | if (calcfast == 1 ) {
2801		Lveil_cie =0 ;
2802		}
2803		Lveil_cie_sum =  Lveil_cie_sum + Lveil_cie;
2804	<	printf("%i %f %f %f %f %.10f %f %f %f %f %f %f %f %f %f %f %f %f %f\n",
2804	>	if (pict_get_z1_gsn(p,x)<0) {
2805	>	act_gsn=(int)(-pict_get_z1_gsn(p,x));
2806	>	}else{
2807	>	act_gsn=0;
2808	>	}
2809	>	printf("%i %f %f %f %f %.10f %f %f %f %f %f %f %f %f %f %f %f %f %f %i \n",
2810		i, pict_get_npix(p, x), pict_get_av_posx(p, x),
2811		pict_get_ysize(p) - pict_get_av_posy(p, x),
2812		pict_get_av_lum(p, x), pict_get_av_omega(p, x),
#	Line 2126 | Line 2814 | if (calcfast == 1 ) {
2814		pict_get_av_posy(p, x),
2815		posindex_2), lum_backg, lum_task,
2816		E_v, E_v_dir, abs_max, sigma * 180 / 3.1415927
2817	<	,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
2130	<
2131	<	);
2817	>	,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 );
2818		}
2819		}
2820		}
#	Line 2168 | Line 2854 | if (calcfast == 1 ) {
2854		}
2855
2856		if (calcfast == 0) {
2857	+	lum_a= E_v2/3.1415927;
2858		dgi = get_dgi(p, lum_backg, igs, posindex_2);
2859		ugr = get_ugr(p, lum_backg, igs, posindex_2);
2860	+	ugp = get_ugp (p,lum_backg, igs, posindex_2);
2861	+	ugr_exp = get_ugr_exp (p,lum_backg_cos,lum_a, igs, posindex_2);
2862	+	dgi_mod = get_dgi_mod(p, lum_a, igs, posindex_2);
2863		cgi = get_cgi(p, E_v, E_v_dir, igs, posindex_2);
2864	<	vcp = get_vcp(p, avlum, igs, posindex_2);
2864	>	dgr = get_dgr(p, avlum, igs, posindex_2);
2865	>	vcp = get_vcp(dgr);
2866		Lveil = get_disability(p, avlum, igs);
2867	+	if (no_glaresources == 0) {
2868	+	dgi = 0.0;
2869	+	ugr = 0.0;
2870	+	ugp = 0.0;
2871	+	ugr_exp = 0.0;
2872	+	dgi_mod = 0.0;
2873	+	cgi = 0.0;
2874	+	dgr = 0.0;
2875	+	vcp = 100.0;
2876	+	}
2877		}
2878		/* check dgp range */
2879		if (dgp <= 0.2) {
#	Line 2198 | Line 2899 | if (calcfast == 0) {
2899		dgp =low_light_corr*dgp;
2900		dgp =age_corr_factor*dgp;
2901		}
2902	+	muc_rvar_get_median(s_noposweight,lum_nopos_median);
2903	+	muc_rvar_get_median(s_posweight,lum_pos_median);
2904	+	muc_rvar_get_median(s_posweight2,lum_pos2_median);
2905	+
2906
2907	+
2908	+
2909		printf
2910	<	("dgp,av_lum,E_v,lum_backg,E_v_dir,dgi,ugr,vcp,cgi,lum_sources,omega_sources,Lveil,Lveil_cie,band_avlum: %f %f %f %f %f %f %f %f %f %f %f %f %f %f \n",
2910	>	("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",
2911		dgp, avlum, E_v, lum_backg, E_v_dir, dgi, ugr, vcp, cgi,
2912	<	lum_sources / omega_sources, omega_sources, Lveil,Lveil_cie_sum,band_avlum);
2912	>	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]);
2913		} else {
2914		if (detail_out2 == 1) {
2915
#	Line 2233 | 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
#	Line 2265 | Line 2979 | has to be re-written from scratch....
2979
2980		/*output  */
2981		if (checkfile == 1) {
2982	+
2983	+
2984		pict_write(p, file_out);
2985		}
2986

Diff Legend

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