--- ray/src/px/pcond4.c	1996/10/11 10:23:59	3.10
+++ ray/src/px/pcond4.c	1997/01/28 16:31:17	3.14
@@ -1,4 +1,4 @@
-/* Copyright (c) 1996 Regents of the University of California */
+/* Copyright (c) 1997 Regents of the University of California */
 
 #ifndef lint
 static char SCCSid[] = "$SunId$ LBL";
@@ -14,13 +14,8 @@ static char SCCSid[] = "$SunId$ LBL";
 
 #define	VADAPT		0.08		/* fraction of adaptation from veil */
 
-extern COLOR	*fovimg;		/* foveal (1 degree) averaged image */
-extern short	fvxr, fvyr;		/* foveal image resolution */
+static COLOR	*veilimg = NULL;	/* veiling image */
 
-#define fovscan(y)	(fovimg+(y)*fvxr)
-
-static COLOR	*veilimg;		/* veiling image */
-
 #define veilscan(y)	(veilimg+(y)*fvxr)
 
 static float	(*raydir)[3] = NULL;	/* ray direction for each pixel */
@@ -83,6 +78,9 @@ compveil()				/* compute veiling image */
 	COLOR	ctmp, vsum;
 	int	px, py;
 	register int	x, y;
+
+	if (veilimg != NULL)		/* already done? */
+		return;
 					/* compute ray directions */
 	compraydir();
 					/* compute veil image */
@@ -100,10 +98,10 @@ compveil()				/* compute veiling image */
 							rdirscan(y)[x]);
 					if (t2 <= FTINY) continue;
 					/*	use approximation instead
-					t2 = acos(t2);
-					t2 = 1./(t2*t2);
+					t3 = acos(t2);
+					t2 = t2/(t3*t3);
 					*/
-					t2 = .5 / (1. - t2);
+					t2 *= .5 / (1. - t2);
 					copycolor(ctmp, fovscan(y)[x]);
 					scalecolor(ctmp, t2);
 					addcolor(vsum, ctmp);
@@ -113,6 +111,13 @@ compveil()				/* compute veiling image */
 			scalecolor(vsum, VADAPT/t2sum);
 			copycolor(veilscan(py)[px], vsum);
 		}
+					/* modify FOV sample image */
+	for (y = 0; y < fvyr; y++)
+		for (x = 0; x < fvxr; x++) {
+			scalecolor(fovscan(y)[x], 1.-VADAPT);
+			addcolor(fovscan(y)[x], veilscan(y)[x]);
+		}
+	comphist();			/* recompute histogram */
 }
 
 
@@ -166,32 +171,9 @@ float	*inpacuD;		/* input acuity data (cycles/degree) 
 double
 hacuity(La)			/* return visual acuity in cycles/degree */
 double	La;
-{				/* data due to S. Shaler (we should fit it!) */
-#define NPOINTS	20
-	static float	l10lum[NPOINTS] = {
-		-3.10503,-2.66403,-2.37703,-2.09303,-1.64403,-1.35803,
-		-1.07403,-0.67203,-0.38503,-0.10103,0.29397,0.58097,0.86497,
-		1.25697,1.54397,1.82797,2.27597,2.56297,2.84697,3.24897
-	};
-	static float	resfreq[NPOINTS] = {
-		2.09,3.28,3.79,4.39,6.11,8.83,10.94,18.66,23.88,31.05,37.42,
-		37.68,41.60,43.16,45.30,47.00,48.43,48.32,51.06,51.09
-	};
-	double	l10La;
-	register int	i;
-					/* check limits */
-	if (La <= 7.85e-4)
-		return(resfreq[0]);
-	if (La >= 1.78e3)
-		return(resfreq[NPOINTS-1]);
-					/* interpolate data */
-	l10La = log10(La);
-	for (i = 0; i < NPOINTS-2 && l10lum[i+1] <= l10La; i++)
-		;
-	return( ( (l10lum[i+1] - l10La)*resfreq[i] +
-			(l10La - l10lum[i])*resfreq[i+1] ) /
-			(l10lum[i+1] - l10lum[i]) );
-#undef NPOINTS
+{
+					/* functional fit */
+	return(17.25*atan(1.4*log10(La) + 0.35) + 25.72);
 }
 
 
@@ -303,11 +285,8 @@ SCANBAR	*sb;
 					/* get scanlines */
 	sl0 = getascan(sb, iy);
 #ifdef DEBUG
-	if (sl0 == NULL) {
-		fprintf(stderr, "%s: internal - cannot backspace in ascanval\n",
-				progname);
-		abort();
-	}
+	if (sl0 == NULL)
+		error(INTERNAL, "cannot backspace in ascanval");
 #endif
 	sl1 = getascan(sb, iy+1);
 					/* 2D linear interpolation */