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

Comparing ray/src/px/pcond4.c (file contents):
Revision 3.5 by greg, Fri Oct 4 18:11:52 1996 UTC vs.
Revision 3.15 by gregl, Wed Aug 13 11:24:48 1997 UTC

# Line 1 | Line 1
1 < /* Copyright (c) 1996 Regents of the University of California */
1 > /* Copyright (c) 1997 Regents of the University of California */
2  
3   #ifndef lint
4   static char SCCSid[] = "$SunId$ LBL";
# Line 14 | Line 14 | static char SCCSid[] = "$SunId$ LBL";
14  
15   #define VADAPT          0.08            /* fraction of adaptation from veil */
16  
17 < extern COLOR    *fovimg;                /* foveal (1 degree) averaged image */
18 < extern short    fvxr, fvyr;             /* foveal image resolution */
17 > static COLOR    *veilimg = NULL;        /* veiling image */
18  
20 #define fovscan(y)      (fovimg+(y)*fvxr)
21
22 static COLOR    *veilimg;               /* veiling image */
23
19   #define veilscan(y)     (veilimg+(y)*fvxr)
20  
21   static float    (*raydir)[3] = NULL;    /* ray direction for each pixel */
# Line 83 | Line 78 | compveil()                             /* compute veiling image */
78          COLOR   ctmp, vsum;
79          int     px, py;
80          register int    x, y;
81 +
82 +        if (veilimg != NULL)            /* already done? */
83 +                return;
84                                          /* compute ray directions */
85          compraydir();
86                                          /* compute veil image */
# Line 100 | Line 98 | compveil()                             /* compute veiling image */
98                                                          rdirscan(y)[x]);
99                                          if (t2 <= FTINY) continue;
100                                          /*      use approximation instead
101 <                                        t2 = acos(t2);
102 <                                        t2 = 1./(t2*t2);
101 >                                        t3 = acos(t2);
102 >                                        t2 = t2/(t3*t3);
103                                          */
104 <                                        t2 = .5 / (1. - t2);
104 >                                        t2 *= .5 / (1. - t2);
105                                          copycolor(ctmp, fovscan(y)[x]);
106                                          scalecolor(ctmp, t2);
107                                          addcolor(vsum, ctmp);
108                                          t2sum += t2;
109                                  }
110                          /* VADAPT of original is subtracted in addveil() */
111 <                        scalecolor(vsum, VADAPT/t2sum);
111 >                        if (t2sum > FTINY)
112 >                                scalecolor(vsum, VADAPT/t2sum);
113                          copycolor(veilscan(py)[px], vsum);
114                  }
115 +                                        /* modify FOV sample image */
116 +        for (y = 0; y < fvyr; y++)
117 +                for (x = 0; x < fvxr; x++) {
118 +                        scalecolor(fovscan(y)[x], 1.-VADAPT);
119 +                        addcolor(fovscan(y)[x], veilscan(y)[x]);
120 +                }
121 +        comphist();                     /* recompute histogram */
122   }
123  
124  
# Line 126 | Line 132 | int    y;
132          register int    x, i;
133  
134          vy = dy = (y+.5)/numscans(&inpres)*fvyr - .5;
135 <        if (vy >= fvyr-1) vy--;
135 >        while (vy >= fvyr-1) vy--;
136          dy -= (double)vy;
137          for (x = 0; x < scanlen(&inpres); x++) {
138                  vx = dx = (x+.5)/scanlen(&inpres)*fvxr - .5;
139 <                if (vx >= fvxr-1) vx--;
139 >                while (vx >= fvxr-1) vx--;
140                  dx -= (double)vx;
141                  for (i = 0; i < 3; i++) {
142                          lv = (1.-dy)*colval(veilscan(vy)[vx],i) +
# Line 146 | Line 152 | int    y;
152  
153   /****************** ACUITY STUFF *******************/
154  
155 < typedef struct scanbar {
156 <        short   sampr;          /* sample area size (power of 2) */
155 > typedef struct {
156 >        short   sampe;          /* sample area size (exponent of 2) */
157          short   nscans;         /* number of scanlines in this bar */
158          int     len;            /* individual scanline length */
153        struct scanbar  *next;  /* next higher resolution scanbar */
159          int     nread;          /* number of scanlines loaded */
160 <                        /* followed by the scanline data */
160 >        COLOR   *sdata;         /* scanbar data */
161   } SCANBAR;
162  
163 < #define bscan(sb,y)     ((COLOR *)((sb)+1)+((y)%(sb)->nscans)*(sb)->len)
163 > #define bscan(sb,y)     ((COLOR *)(sb)->sdata+((y)%(sb)->nscans)*(sb)->len)
164  
165   SCANBAR *rootbar;               /* root scan bar (lowest resolution) */
166  
# Line 167 | Line 172 | float  *inpacuD;               /* input acuity data (cycles/degree)
172   double
173   hacuity(La)                     /* return visual acuity in cycles/degree */
174   double  La;
175 < {                               /* data due to S. Shaler (we should fit it!) */
176 < #define NPOINTS 20
177 <        static float    l10lum[NPOINTS] = {
173 <                -3.10503,-2.66403,-2.37703,-2.09303,-1.64403,-1.35803,
174 <                -1.07403,-0.67203,-0.38503,-0.10103,0.29397,0.58097,0.86497,
175 <                1.25697,1.54397,1.82797,2.27597,2.56297,2.84697,3.24897
176 <        };
177 <        static float    resfreq[NPOINTS] = {
178 <                2.09,3.28,3.79,4.39,6.11,8.83,10.94,18.66,23.88,31.05,37.42,
179 <                37.68,41.60,43.16,45.30,47.00,48.43,48.32,51.06,51.09
180 <        };
181 <        double  l10La;
182 <        register int    i;
183 <                                        /* check limits */
184 <        if (La <= 7.85e-4)
185 <                return(resfreq[0]);
186 <        if (La >= 1.78e3)
187 <                return(resfreq[NPOINTS-1]);
188 <                                        /* interpolate data */
189 <        l10La = log10(La);
190 <        for (i = 0; i < NPOINTS-2 && l10lum[i+1] <= l10La; i++)
191 <                ;
192 <        return( ( (l10lum[i+1] - l10La)*resfreq[i] +
193 <                        (l10La - l10lum[i])*resfreq[i+1] ) /
194 <                        (l10lum[i+1] - l10lum[i]) );
195 < #undef NPOINTS
175 > {
176 >                                        /* functional fit */
177 >        return(17.25*atan(1.4*log10(La) + 0.35) + 25.72);
178   }
179  
180  
# Line 204 | Line 186 | int    y;
186          register COLOR  *sl0, *sl1, *mysl;
187          register int    i;
188  
189 <        if (y < sb->nread - sb->nscans) {
190 <                fprintf(stderr, "%s: internal - cannot backspace in getascan\n",
209 <                                progname);
210 <                exit(1);
211 <        }
189 >        if (y < sb->nread - sb->nscans)                 /* too far back? */
190 >                return(NULL);
191          for ( ; y >= sb->nread; sb->nread++) {          /* read as necessary */
192                  mysl = bscan(sb, sb->nread);
193 <                if (sb->sampr == 1) {
193 >                if (sb->sampe == 0) {
194                          if (freadscan(mysl, sb->len, infp) < 0) {
195                                  fprintf(stderr, "%s: %s: scanline read error\n",
196                                                  progname, infn);
197                                  exit(1);
198                          }
199                  } else {
200 <                        sl0 = getascan(sb->next, 2*y);
201 <                        sl1 = getascan(sb->next, 2*y+1);
200 >                        sl0 = getascan(sb+1, 2*y);
201 >                        if (sl0 == NULL)
202 >                                return(NULL);
203 >                        sl1 = getascan(sb+1, 2*y+1);
204                          for (i = 0; i < sb->len; i++) {
205                                  copycolor(mysl[i], sl0[2*i]);
206                                  addcolor(mysl[i], sl0[2*i+1]);
# Line 243 | Line 224 | int    y;
224          register int    x;
225                                          /* compute foveal y position */
226          iy = dy = (y+.5)/numscans(&inpres)*fvyr - .5;
227 <        if (iy >= fvyr-1) iy--;
227 >        while (iy >= fvyr-1) iy--;
228          dy -= (double)iy;
229          for (x = 0; x < scanlen(&inpres); x++) {
230                                          /* compute foveal x position */
231                  ix = dx = (x+.5)/scanlen(&inpres)*fvxr - .5;
232 <                if (ix >= fvxr-1) ix--;
232 >                while (ix >= fvxr-1) ix--;
233                  dx -= (double)ix;
234                                          /* interpolate sample rate */
235                  sr = (1.-dy)*((1.-dx)*tsampr(ix,iy) + dx*tsampr(ix+1,iy)) +
# Line 268 | Line 249 | double sr;
249          double  d;
250          register SCANBAR        *sb0;
251  
252 <        for (sb0 = rootbar; sb0->next != NULL && sb0->next->sampr > sr;
272 <                        sb0 = sb0->next)
252 >        for (sb0 = rootbar; sb0->sampe != 0 && 1<<sb0[1].sampe > sr; sb0++)
253                  ;
254          ascanval(col, x, y, sb0);
255 <        if (sb0->next == NULL)          /* don't extrapolate highest */
255 >        if (sb0->sampe == 0)            /* don't extrapolate highest */
256                  return;
257 <        ascanval(c1, x, y, sb0->next);
258 <        d = (sb0->sampr - sr)/(sb0->sampr - sb0->next->sampr);
257 >        ascanval(c1, x, y, sb0+1);
258 >        d = ((1<<sb0->sampe) - sr)/(1<<sb0[1].sampe);
259          scalecolor(col, 1.-d);
260          scalecolor(c1, d);
261          addcolor(col, c1);
# Line 291 | Line 271 | SCANBAR        *sb;
271          double  dx, dy;
272          int     ix, iy;
273  
274 <        ix = dx = (x+.5)/sb->sampr - .5;
275 <        if (ix >= sb->len-1) ix--;
274 >        if (sb->sampe == 0) {           /* no need to interpolate */
275 >                sl0 = getascan(sb, y);
276 >                copycolor(col, sl0[x]);
277 >                return;
278 >        }
279 >                                        /* compute coordinates for sb */
280 >        ix = dx = (x+.5)/(1<<sb->sampe) - .5;
281 >        while (ix >= sb->len-1) ix--;
282          dx -= (double)ix;
283 <        iy = dy = (y+.5)/sb->sampr - .5;
284 <        if (iy >= numscans(&inpres)/sb->sampr-1) iy--;
283 >        iy = dy = (y+.5)/(1<<sb->sampe) - .5;
284 >        while (iy >= (numscans(&inpres)>>sb->sampe)-1) iy--;
285          dy -= (double)iy;
286                                          /* get scanlines */
287          sl0 = getascan(sb, iy);
288 + #ifdef DEBUG
289 +        if (sl0 == NULL)
290 +                error(INTERNAL, "cannot backspace in ascanval");
291 + #endif
292          sl1 = getascan(sb, iy+1);
293                                          /* 2D linear interpolation */
294          copycolor(col, sl0[ix]);
# Line 314 | Line 304 | SCANBAR        *sb;
304          scalecolor(col, 1.-dy);
305          scalecolor(c1y, dy);
306          addcolor(col, c1y);
307 +        for (ix = 0; ix < 3; ix++)      /* make sure no negative */
308 +                if (colval(col,ix) < 0.)
309 +                        colval(col,ix) = 0.;
310   }
311  
312  
313   SCANBAR *
314 < sballoc(sr, ns, sl)             /* allocate scanbar */
315 < int     sr;             /* sampling rate */
314 > sballoc(se, ns, sl)             /* allocate scanbar */
315 > int     se;             /* sampling rate exponent */
316   int     ns;             /* number of scanlines */
317   int     sl;             /* original scanline length */
318   {
319 +        SCANBAR *sbarr;
320          register SCANBAR        *sb;
321  
322 <        sb = (SCANBAR *)malloc(sizeof(SCANBAR)+(sl/sr)*ns*sizeof(COLOR));
322 >        sbarr = sb = (SCANBAR *)malloc((se+1)*sizeof(SCANBAR));
323          if (sb == NULL)
324                  syserror("malloc");
325 <        sb->nscans = ns;
326 <        sb->len = sl/sr;
327 <        sb->nread = 0;
328 <        if ((sb->sampr = sr) > 1)
329 <                sb->next = sballoc(sr/2, ns*2, sl);
330 <        else
331 <                sb->next = NULL;
332 <        return(sb);
325 >        do {
326 >                sb->sampe = se;
327 >                sb->len = sl>>se;
328 >                sb->nscans = ns;
329 >                sb->sdata = (COLOR *)malloc(sb->len*ns*sizeof(COLOR));
330 >                if (sb->sdata == NULL)
331 >                        syserror("malloc");
332 >                sb->nread = 0;
333 >                ns <<= 1;
334 >                sb++;
335 >        } while (--se >= 0);
336 >        return(sbarr);
337   }
338  
339  
# Line 363 | Line 361 | initacuity()                   /* initialize variable acuity sampling
361                          }
362                          fcross(cp, diffx, diffy);
363                          omega = 0.5 * sqrt(DOT(cp,cp));
364 <                        if (omega <= FTINY)
364 >                        if (omega <= FTINY*FTINY)
365                                  tsampr(x,y) = 1.;
366                          else if ((tsampr(x,y) = PI/180. / sqrt(omega) /
367                                          hacuity(plum(fovscan(y)[x]))) > maxsr)
# Line 375 | Line 373 | initacuity()                   /* initialize variable acuity sampling
373                  tsampr(x,fvyr-1) = tsampr(x,fvyr-2);
374          }
375          for (y = 0; y < fvyr; y++) {
376 <                tsampr(y,0) = tsampr(y,1);
377 <                tsampr(y,fvxr-1) = tsampr(y,fvxr-2);
376 >                tsampr(0,y) = tsampr(1,y);
377 >                tsampr(fvxr-1,y) = tsampr(fvxr-2,y);
378          }
379                                          /* initialize with next power of two */
380 <        rootbar = sballoc(2<<(int)(log(maxsr)/log(2.)), 2, scanlen(&inpres));
380 >        rootbar = sballoc((int)(log(maxsr)/log(2.))+1, 2, scanlen(&inpres));
381   }

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines