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

Comparing ray/src/rt/ambient.c (file contents):
Revision 1.14 by greg, Fri Jun 7 10:04:29 1991 UTC vs.
Revision 2.92 by greg, Fri Nov 21 00:30:11 2014 UTC

# Line 1 | Line 1
1 /* Copyright (c) 1991 Regents of the University of California */
2
1   #ifndef lint
2 < static char SCCSid[] = "$SunId$ LBL";
2 > static const char       RCSid[] = "$Id$";
3   #endif
6
4   /*
5   *  ambient.c - routines dealing with ambient (inter-reflected) component.
6   *
7 < *  The macro AMBFLUSH (if defined) is the number of ambient values
11 < *      to wait before flushing to the ambient file.
12 < *
13 < *     5/9/86
7 > *  Declarations of external symbols in ambient.h
8   */
9  
10 < #include  "ray.h"
10 > #include "copyright.h"
11  
12 < #include  "octree.h"
12 > #include <string.h>
13  
14 + #include  "platform.h"
15 + #include  "ray.h"
16   #include  "otypes.h"
17 <
17 > #include  "resolu.h"
18   #include  "ambient.h"
23
19   #include  "random.h"
20  
21 < #define  OCTSCALE       0.5     /* ceil((valid rad.)/(cube size)) */
21 > #ifndef  OCTSCALE
22 > #define  OCTSCALE       1.0     /* ceil((valid rad.)/(cube size)) */
23 > #endif
24  
25 < typedef struct ambtree {
29 <        AMBVAL  *alist;         /* ambient value list */
30 <        struct ambtree  *kid;   /* 8 child nodes */
31 < }  AMBTREE;                     /* ambient octree */
25 > extern char  *shm_boundary;     /* memory sharing boundary */
26  
27 < extern CUBE  thescene;          /* contains space boundaries */
27 > #ifndef  MAXASET
28 > #define  MAXASET        4095    /* maximum number of elements in ambient set */
29 > #endif
30 > OBJECT  ambset[MAXASET+1]={0};  /* ambient include/exclude set */
31  
32 < #define  MAXASET        511     /* maximum number of elements in ambient set */
33 < OBJECT  ambset[MAXASET+1]={0};  /* ambient include/exclude set */
32 > double  maxarad;                /* maximum ambient radius */
33 > double  minarad;                /* minimum ambient radius */
34  
35 < double  maxarad;                /* maximum ambient radius */
39 < double  minarad;                /* minimum ambient radius */
35 > static AMBTREE  atrunk;         /* our ambient trunk node */
36  
41 static AMBTREE  atrunk;         /* our ambient trunk node */
42
37   static FILE  *ambfp = NULL;     /* ambient file pointer */
38 + static int  nunflshed = 0;      /* number of unflushed ambient values */
39  
40 < #define  newambval()    (AMBVAL *)bmalloc(sizeof(AMBVAL))
40 > #ifndef SORT_THRESH
41 > #ifdef SMLMEM
42 > #define SORT_THRESH     ((16L<<20)/sizeof(AMBVAL))
43 > #else
44 > #define SORT_THRESH     ((64L<<20)/sizeof(AMBVAL))
45 > #endif
46 > #endif
47 > #ifndef SORT_INTVL
48 > #define SORT_INTVL      (SORT_THRESH<<1)
49 > #endif
50 > #ifndef MAX_SORT_INTVL
51 > #define MAX_SORT_INTVL  (SORT_INTVL<<6)
52 > #endif
53  
47 #define  newambtree()   (AMBTREE *)calloc(8, sizeof(AMBTREE))
54  
55 + static double  avsum = 0.;              /* computed ambient value sum (log) */
56 + static unsigned int  navsum = 0;        /* number of values in avsum */
57 + static unsigned int  nambvals = 0;      /* total number of indirect values */
58 + static unsigned int  nambshare = 0;     /* number of values from file */
59 + static unsigned long  ambclock = 0;     /* ambient access clock */
60 + static unsigned long  lastsort = 0;     /* time of last value sort */
61 + static long  sortintvl = SORT_INTVL;    /* time until next sort */
62 + static FILE  *ambinp = NULL;            /* auxiliary file for input */
63 + static long  lastpos = -1;              /* last flush position */
64  
65 < setambient(afile)                       /* initialize calculation */
66 < char  *afile;
65 > #define MAXACLOCK       (1L<<30)        /* clock turnover value */
66 >        /*
67 >         * Track access times unless we are sharing ambient values
68 >         * through memory on a multiprocessor, when we want to avoid
69 >         * claiming our own memory (copy on write).  Go ahead anyway
70 >         * if more than two thirds of our values are unshared.
71 >         * Compile with -Dtracktime=0 to turn this code off.
72 >         */
73 > #ifndef tracktime
74 > #define tracktime       (shm_boundary == NULL || nambvals > 3*nambshare)
75 > #endif
76 >
77 > #define  AMBFLUSH       (BUFSIZ/AMBVALSIZ)
78 >
79 > #define  newambval()    (AMBVAL *)malloc(sizeof(AMBVAL))
80 > #define  freeav(av)     free((void *)av);
81 >
82 > static void initambfile(int creat);
83 > static void avsave(AMBVAL *av);
84 > static AMBVAL *avstore(AMBVAL  *aval);
85 > static AMBTREE *newambtree(void);
86 > static void freeambtree(AMBTREE  *atp);
87 >
88 > typedef void unloadtf_t(AMBVAL *);
89 > static unloadtf_t avinsert;
90 > static unloadtf_t av2list;
91 > static unloadtf_t avfree;
92 > static void unloadatree(AMBTREE  *at, unloadtf_t *f);
93 >
94 > static int aposcmp(const void *avp1, const void *avp2);
95 > static int avlmemi(AMBVAL *avaddr);
96 > static void sortambvals(int always);
97 >
98 > #ifdef  F_SETLKW
99 > static void aflock(int  typ);
100 > #endif
101 >
102 >
103 > void
104 > setambres(                              /* set ambient resolution */
105 >        int  ar
106 > )
107   {
108 <        long  ftell();
109 <        AMBVAL  amb;
108 >        ambres = ar < 0 ? 0 : ar;               /* may be done already */
109 >                                                /* set min & max radii */
110 >        if (ar <= 0) {
111 >                minarad = 0;
112 >                maxarad = thescene.cusize*0.2;
113 >        } else {
114 >                minarad = thescene.cusize / ar;
115 >                maxarad = 64.0 * minarad;               /* heuristic */
116 >                if (maxarad > thescene.cusize*0.2)
117 >                        maxarad = thescene.cusize*0.2;
118 >        }
119 >        if (minarad <= FTINY)
120 >                minarad = 10.0*FTINY;
121 >        if (maxarad <= minarad)
122 >                maxarad = 64.0 * minarad;
123 > }
124  
56        maxarad = thescene.cusize / 2.0;                /* maximum radius */
57                                                        /* minimum radius */
58        minarad = ambres > 0 ? thescene.cusize/ambres : 0.0;
125  
126 <                                        /* open ambient file */
127 <        if (afile != NULL)
128 <                if ((ambfp = fopen(afile, "r+")) != NULL) {
129 <                        while (fread((char *)&amb,sizeof(AMBVAL),1,ambfp) == 1)
130 <                                avinsert(&amb, &atrunk, thescene.cuorg,
131 <                                                thescene.cusize);
132 <                                                        /* align */
133 <                        fseek(ambfp, -(ftell(ambfp)%sizeof(AMBVAL)), 1);
134 <                } else if ((ambfp = fopen(afile, "w")) == NULL) {
135 <                        sprintf(errmsg, "cannot open ambient file \"%s\"",
136 <                                        afile);
137 <                        error(SYSTEM, errmsg);
126 > void
127 > setambacc(                              /* set ambient accuracy */
128 >        double  newa
129 > )
130 > {
131 >        static double   olda;           /* remember previous setting here */
132 >        
133 >        newa *= (newa > 0);
134 >        if (fabs(newa - olda) >= .05*(newa + olda)) {
135 >                ambacc = newa;
136 >                if (nambvals > 0)
137 >                        sortambvals(1);         /* rebuild tree */
138 >        }
139 > }
140 >
141 >
142 > void
143 > setambient(void)                                /* initialize calculation */
144 > {
145 >        int     readonly = 0;
146 >        long    flen;
147 >        AMBVAL  amb;
148 >                                                /* make sure we're fresh */
149 >        ambdone();
150 >                                                /* init ambient limits */
151 >        setambres(ambres);
152 >        setambacc(ambacc);
153 >        if (ambfile == NULL || !ambfile[0])
154 >                return;
155 >        if (ambacc <= FTINY) {
156 >                sprintf(errmsg, "zero ambient accuracy so \"%s\" not opened",
157 >                                ambfile);
158 >                error(WARNING, errmsg);
159 >                return;
160 >        }
161 >                                                /* open ambient file */
162 >        if ((ambfp = fopen(ambfile, "r+")) == NULL)
163 >                readonly = (ambfp = fopen(ambfile, "r")) != NULL;
164 >        if (ambfp != NULL) {
165 >                initambfile(0);                 /* file exists */
166 >                lastpos = ftell(ambfp);
167 >                while (readambval(&amb, ambfp))
168 >                        avstore(&amb);
169 >                nambshare = nambvals;           /* share loaded values */
170 >                if (readonly) {
171 >                        sprintf(errmsg,
172 >                                "loaded %u values from read-only ambient file",
173 >                                        nambvals);
174 >                        error(WARNING, errmsg);
175 >                        fclose(ambfp);          /* close file so no writes */
176 >                        ambfp = NULL;
177 >                        return;                 /* avoid ambsync() */
178                  }
179 +                                                /* align file pointer */
180 +                lastpos += (long)nambvals*AMBVALSIZ;
181 +                flen = lseek(fileno(ambfp), (off_t)0, SEEK_END);
182 +                if (flen != lastpos) {
183 +                        sprintf(errmsg,
184 +                        "ignoring last %ld values in ambient file (corrupted)",
185 +                                        (flen - lastpos)/AMBVALSIZ);
186 +                        error(WARNING, errmsg);
187 +                        fseek(ambfp, lastpos, SEEK_SET);
188 + #ifndef _WIN32 /* XXX we need a replacement for that one */
189 +                        ftruncate(fileno(ambfp), (off_t)lastpos);
190 + #endif
191 +                }
192 +        } else if ((ambfp = fopen(ambfile, "w+")) != NULL) {
193 +                initambfile(1);                 /* else create new file */
194 +                fflush(ambfp);
195 +                lastpos = ftell(ambfp);
196 +        } else {
197 +                sprintf(errmsg, "cannot open ambient file \"%s\"", ambfile);
198 +                error(SYSTEM, errmsg);
199 +        }
200 + #ifdef getc_unlocked
201 +        flockfile(ambfp);                       /* application-level lock */
202 + #endif
203 + #ifdef  F_SETLKW
204 +        aflock(F_UNLCK);                        /* release file */
205 + #endif
206   }
207  
208  
209 < ambnotify(obj)                  /* record new modifier */
210 < OBJECT  obj;
209 > void
210 > ambdone(void)                   /* close ambient file and free memory */
211   {
212 +        if (ambfp != NULL) {            /* close ambient file */
213 +                ambsync();
214 +                fclose(ambfp);
215 +                ambfp = NULL;
216 +                if (ambinp != NULL) {  
217 +                        fclose(ambinp);
218 +                        ambinp = NULL;
219 +                }
220 +                lastpos = -1;
221 +        }
222 +                                        /* free ambient tree */
223 +        unloadatree(&atrunk, &avfree);
224 +                                        /* reset state variables */
225 +        avsum = 0.;
226 +        navsum = 0;
227 +        nambvals = 0;
228 +        nambshare = 0;
229 +        ambclock = 0;
230 +        lastsort = 0;
231 +        sortintvl = SORT_INTVL;
232 + }
233 +
234 +
235 + void
236 + ambnotify(                      /* record new modifier */
237 +        OBJECT  obj
238 + )
239 + {
240          static int  hitlimit = 0;
241 <        register OBJREC  *o = objptr(obj);
242 <        register char  **amblp;
241 >        OBJREC   *o;
242 >        char  **amblp;
243  
244 +        if (obj == OVOID) {             /* starting over */
245 +                ambset[0] = 0;
246 +                hitlimit = 0;
247 +                return;
248 +        }
249 +        o = objptr(obj);
250          if (hitlimit || !ismodifier(o->otype))
251                  return;
252          for (amblp = amblist; *amblp != NULL; amblp++)
# Line 94 | Line 261 | OBJECT  obj;
261                  }
262   }
263  
264 + /************ THE FOLLOWING ROUTINES DIFFER BETWEEN NEW & OLD ***************/
265  
266 < ambient(acol, r)                /* compute ambient component for ray */
267 < COLOR  acol;
268 < register RAY  *r;
266 > #ifndef OLDAMB
267 >
268 > #define tfunc(lwr, x, upr)      (((x)-(lwr))/((upr)-(lwr)))
269 >
270 > static int      plugaleak(RAY *r, AMBVAL *ap, FVECT anorm, double ang);
271 > static double   sumambient(COLOR acol, RAY *r, FVECT rn, int al,
272 >                                AMBTREE *at, FVECT c0, double s);
273 > static int      makeambient(COLOR acol, RAY *r, FVECT rn, int al);
274 > static int      extambient(COLOR cr, AMBVAL *ap, FVECT pv, FVECT nv,
275 >                                FVECT uvw[3]);
276 >
277 > void
278 > multambient(            /* compute ambient component & multiply by coef. */
279 >        COLOR  aval,
280 >        RAY  *r,
281 >        FVECT  nrm
282 > )
283   {
284          static int  rdepth = 0;                 /* ambient recursion */
285 <        double  wsum;
285 >        COLOR   acol;
286 >        int     ok;
287 >        double  d, l;
288  
105        rdepth++;                               /* increment level */
106
289          if (ambdiv <= 0)                        /* no ambient calculation */
290                  goto dumbamb;
291                                                  /* check number of bounces */
292 <        if (rdepth > ambounce)
292 >        if (rdepth >= ambounce)
293                  goto dumbamb;
294                                                  /* check ambient list */
295          if (ambincl != -1 && r->ro != NULL &&
# Line 115 | Line 297 | register RAY  *r;
297                  goto dumbamb;
298  
299          if (ambacc <= FTINY) {                  /* no ambient storage */
300 <                if (doambient(acol, r, NULL, NULL) == 0.0)
300 >                copycolor(acol, aval);
301 >                rdepth++;
302 >                ok = doambient(acol, r, r->rweight,
303 >                                NULL, NULL, NULL, NULL, NULL);
304 >                rdepth--;
305 >                if (!ok)
306                          goto dumbamb;
307 <                goto done;
307 >                copycolor(aval, acol);
308 >                return;
309          }
310 <                                                /* get ambient value */
310 >
311 >        if (tracktime)                          /* sort to minimize thrashing */
312 >                sortambvals(0);
313 >                                                /* interpolate ambient value */
314          setcolor(acol, 0.0, 0.0, 0.0);
315 <        wsum = sumambient(acol, r, &atrunk, thescene.cuorg, thescene.cusize);
316 <        if (wsum > FTINY)
317 <                scalecolor(acol, 1.0/wsum);
318 <        else if (makeambient(acol, r) == 0.0)
315 >        d = sumambient(acol, r, nrm, rdepth,
316 >                        &atrunk, thescene.cuorg, thescene.cusize);
317 >        if (d > FTINY) {
318 >                d = 1.0/d;
319 >                scalecolor(acol, d);
320 >                multcolor(aval, acol);
321 >                return;
322 >        }
323 >        rdepth++;                               /* need to cache new value */
324 >        ok = makeambient(acol, r, nrm, rdepth-1);
325 >        rdepth--;
326 >        if (ok) {
327 >                multcolor(aval, acol);          /* computed new value */
328 >                return;
329 >        }
330 > dumbamb:                                        /* return global value */
331 >        if ((ambvwt <= 0) | (navsum == 0)) {
332 >                multcolor(aval, ambval);
333 >                return;
334 >        }
335 >        l = bright(ambval);                     /* average in computations */
336 >        if (l > FTINY) {
337 >                d = (log(l)*(double)ambvwt + avsum) /
338 >                                (double)(ambvwt + navsum);
339 >                d = exp(d) / l;
340 >                scalecolor(aval, d);
341 >                multcolor(aval, ambval);        /* apply color of ambval */
342 >        } else {
343 >                d = exp( avsum / (double)navsum );
344 >                scalecolor(aval, d);            /* neutral color */
345 >        }
346 > }
347 >
348 >
349 > /* Plug a potential leak where ambient cache value is occluded */
350 > static int
351 > plugaleak(RAY *r, AMBVAL *ap, FVECT anorm, double ang)
352 > {
353 >        const double    cost70sq = 0.1169778;   /* cos(70deg)^2 */
354 >        RAY             rtst;
355 >        FVECT           vdif;
356 >        double          normdot, ndotd, nadotd;
357 >        double          a, b, c, t[2];
358 >
359 >        ang += 2.*PI*(ang < 0);                 /* check direction flags */
360 >        if ( !(ap->corral>>(int)(ang*(16./PI)) & 1) )
361 >                return(0);
362 >        /*
363 >         * Generate test ray, targeting 20 degrees above sample point plane
364 >         * along surface normal from cache position.  This should be high
365 >         * enough to miss local geometry we don't really care about.
366 >         */
367 >        VSUB(vdif, ap->pos, r->rop);
368 >        normdot = DOT(anorm, r->ron);
369 >        ndotd = DOT(vdif, r->ron);
370 >        nadotd = DOT(vdif, anorm);
371 >        a = normdot*normdot - cost70sq;
372 >        b = 2.0*(normdot*ndotd - nadotd*cost70sq);
373 >        c = ndotd*ndotd - DOT(vdif,vdif)*cost70sq;
374 >        if (quadratic(t, a, b, c) != 2)
375 >                return(1);                      /* should rarely happen */
376 >        if (t[1] <= FTINY)
377 >                return(0);                      /* should fail behind test */
378 >        rayorigin(&rtst, SHADOW, r, NULL);
379 >        VSUM(rtst.rdir, vdif, anorm, t[1]);     /* further dist. > plane */
380 >        rtst.rmax = normalize(rtst.rdir);       /* short ray test */
381 >        while (localhit(&rtst, &thescene)) {    /* check for occluder */
382 >                if (rtst.ro->omod != OVOID &&
383 >                                (rtst.clipset == NULL ||
384 >                                        !inset(rtst.clipset, rtst.ro->omod)))
385 >                        return(1);              /* plug light leak */
386 >                VCOPY(rtst.rorg, rtst.rop);     /* skip invisible surface */
387 >                rtst.rmax -= rtst.rot;
388 >                rayclear(&rtst);
389 >        }
390 >        return(0);                              /* seems we're OK */
391 > }
392 >
393 >
394 > static double
395 > sumambient(             /* get interpolated ambient value */
396 >        COLOR  acol,
397 >        RAY  *r,
398 >        FVECT  rn,
399 >        int  al,
400 >        AMBTREE  *at,
401 >        FVECT  c0,
402 >        double  s
403 > )
404 > {                       /* initial limit is 10 degrees plus ambacc radians */
405 >        const double    minangle = 10.0 * PI/180.;
406 >        double          maxangle = minangle + ambacc;
407 >        double          wsum = 0.0;
408 >        FVECT           ck0;
409 >        int             i, j;
410 >        AMBVAL          *av;
411 >
412 >        if (at->kid != NULL) {          /* sum children first */                                
413 >                s *= 0.5;
414 >                for (i = 0; i < 8; i++) {
415 >                        for (j = 0; j < 3; j++) {
416 >                                ck0[j] = c0[j];
417 >                                if (1<<j & i)
418 >                                        ck0[j] += s;
419 >                                if (r->rop[j] < ck0[j] - OCTSCALE*s)
420 >                                        break;
421 >                                if (r->rop[j] > ck0[j] + (1.0+OCTSCALE)*s)
422 >                                        break;
423 >                        }
424 >                        if (j == 3)
425 >                                wsum += sumambient(acol, r, rn, al,
426 >                                                        at->kid+i, ck0, s);
427 >                }
428 >                                        /* good enough? */
429 >                if (wsum >= 0.05 && s > minarad*10.0)
430 >                        return(wsum);
431 >        }
432 >                                        /* adjust maximum angle */
433 >        if (at->alist != NULL && (at->alist->lvl <= al) & (r->rweight < 0.6))
434 >                maxangle = (maxangle - PI/2.)*pow(r->rweight,0.13) + PI/2.;
435 >                                        /* sum this node */
436 >        for (av = at->alist; av != NULL; av = av->next) {
437 >                double  u, v, d, delta_r2, delta_t2;
438 >                COLOR   ct;
439 >                FVECT   uvw[3];
440 >                                        /* record access */
441 >                if (tracktime)
442 >                        av->latick = ambclock;
443 >                /*
444 >                 *  Ambient level test
445 >                 */
446 >                if (av->lvl > al ||     /* list sorted, so this works */
447 >                                (av->lvl == al) & (av->weight < 0.9*r->rweight))
448 >                        break;
449 >                /*
450 >                 *  Direction test using unperturbed normal
451 >                 */
452 >                decodedir(uvw[2], av->ndir);
453 >                d = DOT(uvw[2], r->ron);
454 >                if (d <= 0.0)           /* >= 90 degrees */
455 >                        continue;
456 >                delta_r2 = 2.0 - 2.0*d; /* approx. radians^2 */
457 >                if (delta_r2 >= maxangle*maxangle)
458 >                        continue;
459 >                /*
460 >                 *  Modified ray behind test
461 >                 */
462 >                VSUB(ck0, r->rop, av->pos);
463 >                d = DOT(ck0, uvw[2]);
464 >                if (d < -minarad*ambacc-.001)
465 >                        continue;
466 >                d /= av->rad[0];
467 >                delta_t2 = d*d;
468 >                if (delta_t2 >= ambacc*ambacc)
469 >                        continue;
470 >                /*
471 >                 *  Elliptical radii test based on Hessian
472 >                 */
473 >                decodedir(uvw[0], av->udir);
474 >                VCROSS(uvw[1], uvw[2], uvw[0]);
475 >                d = (u = DOT(ck0, uvw[0])) / av->rad[0];
476 >                delta_t2 += d*d;
477 >                d = (v = DOT(ck0, uvw[1])) / av->rad[1];
478 >                delta_t2 += d*d;
479 >                if (delta_t2 >= ambacc*ambacc)
480 >                        continue;
481 >                /*
482 >                 *  Test for potential light leak
483 >                 */
484 >                if (av->corral && plugaleak(r, av, uvw[2], atan2a(v,u)))
485 >                        continue;
486 >                /*
487 >                 *  Extrapolate value and compute final weight (hat function)
488 >                 */
489 >                if (!extambient(ct, av, r->rop, rn, uvw))
490 >                        continue;
491 >                d = tfunc(maxangle, sqrt(delta_r2), 0.0) *
492 >                        tfunc(ambacc, sqrt(delta_t2), 0.0);
493 >                scalecolor(ct, d);
494 >                addcolor(acol, ct);
495 >                wsum += d;
496 >        }
497 >        return(wsum);
498 > }
499 >
500 >
501 > static int
502 > makeambient(            /* make a new ambient value for storage */
503 >        COLOR  acol,
504 >        RAY  *r,
505 >        FVECT  rn,
506 >        int  al
507 > )
508 > {
509 >        AMBVAL  amb;
510 >        FVECT   uvw[3];
511 >        int     i;
512 >
513 >        amb.weight = 1.0;                       /* compute weight */
514 >        for (i = al; i-- > 0; )
515 >                amb.weight *= AVGREFL;
516 >        if (r->rweight < 0.1*amb.weight)        /* heuristic override */
517 >                amb.weight = 1.25*r->rweight;
518 >        setcolor(acol, AVGREFL, AVGREFL, AVGREFL);
519 >                                                /* compute ambient */
520 >        i = doambient(acol, r, amb.weight,
521 >                        uvw, amb.rad, amb.gpos, amb.gdir, &amb.corral);
522 >        scalecolor(acol, 1./AVGREFL);           /* undo assumed reflectance */
523 >        if (i <= 0 || amb.rad[0] <= FTINY)      /* no Hessian or zero radius */
524 >                return(i);
525 >                                                /* store value */
526 >        VCOPY(amb.pos, r->rop);
527 >        amb.ndir = encodedir(r->ron);
528 >        amb.udir = encodedir(uvw[0]);
529 >        amb.lvl = al;
530 >        copycolor(amb.val, acol);
531 >                                                /* insert into tree */
532 >        avsave(&amb);                           /* and save to file */
533 >        if (rn != r->ron) {                     /* texture */
534 >                VCOPY(uvw[2], r->ron);
535 >                extambient(acol, &amb, r->rop, rn, uvw);
536 >        }
537 >        return(1);
538 > }
539 >
540 >
541 > static int
542 > extambient(             /* extrapolate value at pv, nv */
543 >        COLOR  cr,
544 >        AMBVAL   *ap,
545 >        FVECT  pv,
546 >        FVECT  nv,
547 >        FVECT  uvw[3]
548 > )
549 > {
550 >        static FVECT    my_uvw[3];
551 >        FVECT           v1;
552 >        int             i;
553 >        double          d = 1.0;        /* zeroeth order */
554 >
555 >        if (uvw == NULL) {              /* need local coordinates? */
556 >                decodedir(my_uvw[2], ap->ndir);
557 >                decodedir(my_uvw[0], ap->udir);
558 >                VCROSS(my_uvw[1], my_uvw[2], my_uvw[0]);
559 >                uvw = my_uvw;
560 >        }
561 >        for (i = 3; i--; )              /* gradient due to translation */
562 >                d += (pv[i] - ap->pos[i]) *
563 >                        (ap->gpos[0]*uvw[0][i] + ap->gpos[1]*uvw[1][i]);
564 >
565 >        VCROSS(v1, uvw[2], nv);         /* gradient due to rotation */
566 >        for (i = 3; i--; )
567 >                d += v1[i] * (ap->gdir[0]*uvw[0][i] + ap->gdir[1]*uvw[1][i]);
568 >        
569 >        if (d <= 0.0) {
570 >                setcolor(cr, 0.0, 0.0, 0.0);
571 >                return(0);              /* should not use if we can avoid it */
572 >        }
573 >        copycolor(cr, ap->val);
574 >        scalecolor(cr, d);
575 >        return(1);
576 > }
577 >
578 >
579 > static void
580 > avinsert(                               /* insert ambient value in our tree */
581 >        AMBVAL *av
582 > )
583 > {
584 >        AMBTREE  *at;
585 >        AMBVAL  *ap;
586 >        AMBVAL  avh;
587 >        FVECT  ck0;
588 >        double  s;
589 >        int  branch;
590 >        int  i;
591 >
592 >        if (av->rad[0] <= FTINY)
593 >                error(CONSISTENCY, "zero ambient radius in avinsert");
594 >        at = &atrunk;
595 >        VCOPY(ck0, thescene.cuorg);
596 >        s = thescene.cusize;
597 >        while (s*(OCTSCALE/2) > av->rad[1]*ambacc) {
598 >                if (at->kid == NULL)
599 >                        if ((at->kid = newambtree()) == NULL)
600 >                                error(SYSTEM, "out of memory in avinsert");
601 >                s *= 0.5;
602 >                branch = 0;
603 >                for (i = 0; i < 3; i++)
604 >                        if (av->pos[i] > ck0[i] + s) {
605 >                                ck0[i] += s;
606 >                                branch |= 1 << i;
607 >                        }
608 >                at = at->kid + branch;
609 >        }
610 >        avh.next = at->alist;           /* order by increasing level */
611 >        for (ap = &avh; ap->next != NULL; ap = ap->next)
612 >                if ( ap->next->lvl > av->lvl ||
613 >                                (ap->next->lvl == av->lvl) &
614 >                                (ap->next->weight <= av->weight) )
615 >                        break;
616 >        av->next = ap->next;
617 >        ap->next = (AMBVAL*)av;
618 >        at->alist = avh.next;
619 > }
620 >
621 >
622 > #else /* ! NEWAMB */
623 >
624 > static double   sumambient(COLOR acol, RAY *r, FVECT rn, int al,
625 >                                AMBTREE *at, FVECT c0, double s);
626 > static double   makeambient(COLOR acol, RAY *r, FVECT rn, int al);
627 > static void     extambient(COLOR cr, AMBVAL *ap, FVECT pv, FVECT nv);
628 >
629 >
630 > void
631 > multambient(            /* compute ambient component & multiply by coef. */
632 >        COLOR  aval,
633 >        RAY  *r,
634 >        FVECT  nrm
635 > )
636 > {
637 >        static int  rdepth = 0;                 /* ambient recursion */
638 >        COLOR   acol;
639 >        double  d, l;
640 >
641 >        if (ambdiv <= 0)                        /* no ambient calculation */
642                  goto dumbamb;
643 <        goto done;
643 >                                                /* check number of bounces */
644 >        if (rdepth >= ambounce)
645 >                goto dumbamb;
646 >                                                /* check ambient list */
647 >        if (ambincl != -1 && r->ro != NULL &&
648 >                        ambincl != inset(ambset, r->ro->omod))
649 >                goto dumbamb;
650  
651 < dumbamb:                                        /* return global value */
652 <        copycolor(acol, ambval);
653 < done:                                           /* must finish here! */
651 >        if (ambacc <= FTINY) {                  /* no ambient storage */
652 >                copycolor(acol, aval);
653 >                rdepth++;
654 >                d = doambient(acol, r, r->rweight, NULL, NULL);
655 >                rdepth--;
656 >                if (d <= FTINY)
657 >                        goto dumbamb;
658 >                copycolor(aval, acol);
659 >                return;
660 >        }
661 >
662 >        if (tracktime)                          /* sort to minimize thrashing */
663 >                sortambvals(0);
664 >                                                /* interpolate ambient value */
665 >        setcolor(acol, 0.0, 0.0, 0.0);
666 >        d = sumambient(acol, r, nrm, rdepth,
667 >                        &atrunk, thescene.cuorg, thescene.cusize);
668 >        if (d > FTINY) {
669 >                d = 1.0/d;
670 >                scalecolor(acol, d);
671 >                multcolor(aval, acol);
672 >                return;
673 >        }
674 >        rdepth++;                               /* need to cache new value */
675 >        d = makeambient(acol, r, nrm, rdepth-1);
676          rdepth--;
677 +        if (d > FTINY) {
678 +                multcolor(aval, acol);          /* got new value */
679 +                return;
680 +        }
681 + dumbamb:                                        /* return global value */
682 +        if ((ambvwt <= 0) | (navsum == 0)) {
683 +                multcolor(aval, ambval);
684 +                return;
685 +        }
686 +        l = bright(ambval);                     /* average in computations */
687 +        if (l > FTINY) {
688 +                d = (log(l)*(double)ambvwt + avsum) /
689 +                                (double)(ambvwt + navsum);
690 +                d = exp(d) / l;
691 +                scalecolor(aval, d);
692 +                multcolor(aval, ambval);        /* apply color of ambval */
693 +        } else {
694 +                d = exp( avsum / (double)navsum );
695 +                scalecolor(aval, d);            /* neutral color */
696 +        }
697   }
698  
699  
700 < double
701 < sumambient(acol, r, at, c0, s)          /* get interpolated ambient value */
702 < COLOR  acol;
703 < register RAY  *r;
704 < AMBTREE  *at;
705 < FVECT  c0;
706 < double  s;
700 > static double
701 > sumambient(     /* get interpolated ambient value */
702 >        COLOR  acol,
703 >        RAY  *r,
704 >        FVECT  rn,
705 >        int  al,
706 >        AMBTREE  *at,
707 >        FVECT  c0,
708 >        double  s
709 > )
710   {
711 <        extern double  sqrt();
147 <        double  d, e1, e2, wt, wsum;
711 >        double  d, e1, e2, wt, wsum;
712          COLOR  ct;
713          FVECT  ck0;
714          int  i;
715 <        register int  j;
716 <        register AMBVAL  *av;
717 <                                        /* do this node */
715 >        int  j;
716 >        AMBVAL   *av;
717 >
718          wsum = 0.0;
719 +                                        /* do this node */
720          for (av = at->alist; av != NULL; av = av->next) {
721 +                double  rn_dot = -2.0;
722 +                if (tracktime)
723 +                        av->latick = ambclock;
724                  /*
725 <                 *  Ray strength test.
725 >                 *  Ambient level test.
726                   */
727 <                if (av->lvl > r->rlvl || av->weight < r->rweight-FTINY)
728 <                        continue;
727 >                if (av->lvl > al ||     /* list sorted, so this works */
728 >                                (av->lvl == al) & (av->weight < 0.9*r->rweight))
729 >                        break;
730                  /*
731                   *  Ambient radius test.
732                   */
733 <                e1 = 0.0;
734 <                for (j = 0; j < 3; j++) {
166 <                        d = av->pos[j] - r->rop[j];
167 <                        e1 += d * d;
168 <                }
169 <                e1 /= av->rad * av->rad;
733 >                VSUB(ck0, av->pos, r->rop);
734 >                e1 = DOT(ck0, ck0) / (av->rad * av->rad);
735                  if (e1 > ambacc*ambacc*1.21)
736                          continue;
737                  /*
738 <                 *  Normal direction test.
738 >                 *  Direction test using closest normal.
739                   */
740 <                e2 = (1.0 - DOT(av->dir, r->ron)) * r->rweight;
741 <                if (e2 < 0.0) e2 = 0.0;
742 <                if (e1 + e2 > ambacc*ambacc*1.21)
740 >                d = DOT(av->dir, r->ron);
741 >                if (rn != r->ron) {
742 >                        rn_dot = DOT(av->dir, rn);
743 >                        if (rn_dot > 1.0-FTINY)
744 >                                rn_dot = 1.0-FTINY;
745 >                        if (rn_dot >= d-FTINY) {
746 >                                d = rn_dot;
747 >                                rn_dot = -2.0;
748 >                        }
749 >                }
750 >                e2 = (1.0 - d) * r->rweight;
751 >                if (e2 < 0.0)
752 >                        e2 = 0.0;
753 >                else if (e1 + e2 > ambacc*ambacc*1.21)
754                          continue;
755                  /*
756                   *  Ray behind test.
# Line 183 | Line 759 | double  s;
759                  for (j = 0; j < 3; j++)
760                          d += (r->rop[j] - av->pos[j]) *
761                                          (av->dir[j] + r->ron[j]);
762 <                if (d*0.5 < -minarad*ambacc)
762 >                if (d*0.5 < -minarad*ambacc-.001)
763                          continue;
764                  /*
765                   *  Jittering final test reduces image artifacts.
766                   */
767 <                wt = sqrt(e1) + sqrt(e2);
768 <                wt *= .9 + .2*frandom();
769 <                if (wt > ambacc)
767 >                e1 = sqrt(e1);
768 >                e2 = sqrt(e2);
769 >                wt = e1 + e2;
770 >                if (wt > ambacc*(.9+.2*urand(9015+samplendx)))
771                          continue;
772 +                /*
773 +                 *  Recompute directional error using perturbed normal
774 +                 */
775 +                if (rn_dot > 0.0) {
776 +                        e2 = sqrt((1.0 - rn_dot)*r->rweight);
777 +                        wt = e1 + e2;
778 +                }
779                  if (wt <= 1e-3)
780                          wt = 1e3;
781                  else
782                          wt = 1.0 / wt;
783                  wsum += wt;
784 <                copycolor(ct, av->val);
784 >                extambient(ct, av, r->rop, rn);
785                  scalecolor(ct, wt);
786                  addcolor(acol, ct);
787          }
# Line 216 | Line 800 | double  s;
800                                  break;
801                  }
802                  if (j == 3)
803 <                        wsum += sumambient(acol, r, at->kid+i, ck0, s);
803 >                        wsum += sumambient(acol, r, rn, al,
804 >                                                at->kid+i, ck0, s);
805          }
806          return(wsum);
807   }
808  
809  
810 < double
811 < makeambient(acol, r)            /* make a new ambient value */
812 < COLOR  acol;
813 < register RAY  *r;
810 > static double
811 > makeambient(            /* make a new ambient value for storage */
812 >        COLOR  acol,
813 >        RAY  *r,
814 >        FVECT  rn,
815 >        int  al
816 > )
817   {
818 <        AMBVAL  amb;
818 >        AMBVAL  amb;
819          FVECT   gp, gd;
820 +        int     i;
821  
822 <        amb.rad = doambient(acol, r, gp, gd);   /* compute ambient */
823 <        if (amb.rad == 0.0)
822 >        amb.weight = 1.0;                       /* compute weight */
823 >        for (i = al; i-- > 0; )
824 >                amb.weight *= AVGREFL;
825 >        if (r->rweight < 0.1*amb.weight)        /* heuristic override */
826 >                amb.weight = 1.25*r->rweight;
827 >        setcolor(acol, AVGREFL, AVGREFL, AVGREFL);
828 >                                                /* compute ambient */
829 >        amb.rad = doambient(acol, r, amb.weight, gp, gd);
830 >        if (amb.rad <= FTINY) {
831 >                setcolor(acol, 0.0, 0.0, 0.0);
832                  return(0.0);
833 <                                                /* store it */
833 >        }
834 >        scalecolor(acol, 1./AVGREFL);           /* undo assumed reflectance */
835 >                                                /* store value */
836          VCOPY(amb.pos, r->rop);
837          VCOPY(amb.dir, r->ron);
838 <        amb.lvl = r->rlvl;
240 <        amb.weight = r->rweight;
838 >        amb.lvl = al;
839          copycolor(amb.val, acol);
840          VCOPY(amb.gpos, gp);
841          VCOPY(amb.gdir, gd);
842                                                  /* insert into tree */
843 <        avinsert(&amb, &atrunk, thescene.cuorg, thescene.cusize);
844 <        avsave(&amb);                           /* write to file */
843 >        avsave(&amb);                           /* and save to file */
844 >        if (rn != r->ron)
845 >                extambient(acol, &amb, r->rop, rn);     /* texture */
846          return(amb.rad);
847   }
848  
849  
850 < static
851 < avsave(av)                              /* save an ambient value */
852 < AMBVAL  *av;
850 > static void
851 > extambient(             /* extrapolate value at pv, nv */
852 >        COLOR  cr,
853 >        AMBVAL   *ap,
854 >        FVECT  pv,
855 >        FVECT  nv
856 > )
857   {
858 < #ifdef  AMBFLUSH
859 <        static int  nunflshed = 0;
860 < #endif
861 <        if (ambfp == NULL)
858 >        FVECT  v1;
859 >        int  i;
860 >        double  d;
861 >
862 >        d = 1.0;                        /* zeroeth order */
863 >                                        /* gradient due to translation */
864 >        for (i = 0; i < 3; i++)
865 >                d += ap->gpos[i]*(pv[i]-ap->pos[i]);
866 >                                        /* gradient due to rotation */
867 >        VCROSS(v1, ap->dir, nv);
868 >        d += DOT(ap->gdir, v1);
869 >        if (d <= 0.0) {
870 >                setcolor(cr, 0.0, 0.0, 0.0);
871                  return;
260        if (fwrite((char *)av, sizeof(AMBVAL), 1, ambfp) != 1)
261                goto writerr;
262 #ifdef  AMBFLUSH
263        if (++nunflshed >= AMBFLUSH) {
264                if (fflush(ambfp) == EOF)
265                        goto writerr;
266                nunflshed = 0;
872          }
873 < #endif
874 <        return;
270 < writerr:
271 <        error(SYSTEM, "error writing ambient file");
873 >        copycolor(cr, ap->val);
874 >        scalecolor(cr, d);
875   }
876  
877  
878 < static
879 < avinsert(aval, at, c0, s)               /* insert ambient value in a tree */
880 < AMBVAL  *aval;
881 < register AMBTREE  *at;
279 < FVECT  c0;
280 < double  s;
878 > static void
879 > avinsert(                               /* insert ambient value in our tree */
880 >        AMBVAL *av
881 > )
882   {
883 +        AMBTREE  *at;
884 +        AMBVAL  *ap;
885 +        AMBVAL  avh;
886          FVECT  ck0;
887 +        double  s;
888          int  branch;
889 <        register AMBVAL  *av;
285 <        register int  i;
889 >        int  i;
890  
891 <        if ((av = newambval()) == NULL)
892 <                goto memerr;
893 <        copystruct(av, aval);
894 <        VCOPY(ck0, c0);
891 >        if (av->rad <= FTINY)
892 >                error(CONSISTENCY, "zero ambient radius in avinsert");
893 >        at = &atrunk;
894 >        VCOPY(ck0, thescene.cuorg);
895 >        s = thescene.cusize;
896          while (s*(OCTSCALE/2) > av->rad*ambacc) {
897                  if (at->kid == NULL)
898                          if ((at->kid = newambtree()) == NULL)
899 <                                goto memerr;
899 >                                error(SYSTEM, "out of memory in avinsert");
900                  s *= 0.5;
901                  branch = 0;
902                  for (i = 0; i < 3; i++)
# Line 301 | Line 906 | double  s;
906                          }
907                  at = at->kid + branch;
908          }
909 <        av->next = at->alist;
910 <        at->alist = av;
909 >        avh.next = at->alist;           /* order by increasing level */
910 >        for (ap = &avh; ap->next != NULL; ap = ap->next)
911 >                if ( ap->next->lvl > av->lvl ||
912 >                                (ap->next->lvl == av->lvl) &
913 >                                (ap->next->weight <= av->weight) )
914 >                        break;
915 >        av->next = ap->next;
916 >        ap->next = (AMBVAL*)av;
917 >        at->alist = avh.next;
918 > }
919 >
920 > #endif  /* ! NEWAMB */
921 >
922 > /************* FOLLOWING ROUTINES SAME FOR NEW & OLD METHODS ***************/
923 >
924 > static void
925 > initambfile(            /* initialize ambient file */
926 >        int  cre8
927 > )
928 > {
929 >        extern char  *progname, *octname;
930 >        static char  *mybuf = NULL;
931 >
932 > #ifdef  F_SETLKW
933 >        aflock(cre8 ? F_WRLCK : F_RDLCK);
934 > #endif
935 >        SET_FILE_BINARY(ambfp);
936 >        if (mybuf == NULL)
937 >                mybuf = (char *)bmalloc(BUFSIZ+8);
938 >        setbuf(ambfp, mybuf);
939 >        if (cre8) {                     /* new file */
940 >                newheader("RADIANCE", ambfp);
941 >                fprintf(ambfp, "%s -av %g %g %g -aw %d -ab %d -aa %g ",
942 >                                progname, colval(ambval,RED),
943 >                                colval(ambval,GRN), colval(ambval,BLU),
944 >                                ambvwt, ambounce, ambacc);
945 >                fprintf(ambfp, "-ad %d -as %d -ar %d ",
946 >                                ambdiv, ambssamp, ambres);
947 >                if (octname != NULL)
948 >                        fputs(octname, ambfp);
949 >                fputc('\n', ambfp);
950 >                fprintf(ambfp, "SOFTWARE= %s\n", VersionID);
951 >                fputnow(ambfp);
952 >                fputformat(AMBFMT, ambfp);
953 >                fputc('\n', ambfp);
954 >                putambmagic(ambfp);
955 >        } else if (checkheader(ambfp, AMBFMT, NULL) < 0 || !hasambmagic(ambfp))
956 >                error(USER, "bad ambient file");
957 > }
958 >
959 >
960 > static void
961 > avsave(                         /* insert and save an ambient value */
962 >        AMBVAL  *av
963 > )
964 > {
965 >        avstore(av);
966 >        if (ambfp == NULL)
967 >                return;
968 >        if (writambval(av, ambfp) < 0)
969 >                goto writerr;
970 >        if (++nunflshed >= AMBFLUSH)
971 >                if (ambsync() == EOF)
972 >                        goto writerr;
973          return;
974 < memerr:
975 <        error(SYSTEM, "out of memory in avinsert");
974 > writerr:
975 >        error(SYSTEM, "error writing to ambient file");
976   }
977 +
978 +
979 + static AMBVAL *
980 + avstore(                                /* allocate memory and save aval */
981 +        AMBVAL  *aval
982 + )
983 + {
984 +        AMBVAL  *av;
985 +        double  d;
986 +
987 +        if ((av = newambval()) == NULL)
988 +                error(SYSTEM, "out of memory in avstore");
989 +        *av = *aval;
990 +        av->latick = ambclock;
991 +        av->next = NULL;
992 +        nambvals++;
993 +        d = bright(av->val);
994 +        if (d > FTINY) {                /* add to log sum for averaging */
995 +                avsum += log(d);
996 +                navsum++;
997 +        }
998 +        avinsert(av);                   /* insert in our cache tree */
999 +        return(av);
1000 + }
1001 +
1002 +
1003 + #define ATALLOCSZ       512             /* #/8 trees to allocate at once */
1004 +
1005 + static AMBTREE  *atfreelist = NULL;     /* free ambient tree structures */
1006 +
1007 +
1008 + static AMBTREE *
1009 + newambtree(void)                                /* allocate 8 ambient tree structs */
1010 + {
1011 +        AMBTREE  *atp, *upperlim;
1012 +
1013 +        if (atfreelist == NULL) {       /* get more nodes */
1014 +                atfreelist = (AMBTREE *)malloc(ATALLOCSZ*8*sizeof(AMBTREE));
1015 +                if (atfreelist == NULL)
1016 +                        return(NULL);
1017 +                                        /* link new free list */
1018 +                upperlim = atfreelist + 8*(ATALLOCSZ-1);
1019 +                for (atp = atfreelist; atp < upperlim; atp += 8)
1020 +                        atp->kid = atp + 8;
1021 +                atp->kid = NULL;
1022 +        }
1023 +        atp = atfreelist;
1024 +        atfreelist = atp->kid;
1025 +        memset((char *)atp, '\0', 8*sizeof(AMBTREE));
1026 +        return(atp);
1027 + }
1028 +
1029 +
1030 + static void
1031 + freeambtree(                    /* free 8 ambient tree structs */
1032 +        AMBTREE  *atp
1033 + )
1034 + {
1035 +        atp->kid = atfreelist;
1036 +        atfreelist = atp;
1037 + }
1038 +
1039 +
1040 + static void
1041 + unloadatree(                    /* unload an ambient value tree */
1042 +        AMBTREE  *at,
1043 +        unloadtf_t *f
1044 + )
1045 + {
1046 +        AMBVAL  *av;
1047 +        int  i;
1048 +                                        /* transfer values at this node */
1049 +        for (av = at->alist; av != NULL; av = at->alist) {
1050 +                at->alist = av->next;
1051 +                (*f)(av);
1052 +        }
1053 +        if (at->kid == NULL)
1054 +                return;
1055 +        for (i = 0; i < 8; i++)         /* transfer and free children */
1056 +                unloadatree(at->kid+i, f);
1057 +        freeambtree(at->kid);
1058 +        at->kid = NULL;
1059 + }
1060 +
1061 +
1062 + static struct avl {
1063 +        AMBVAL  *p;
1064 +        unsigned long   t;
1065 + }       *avlist1;                       /* ambient value list with ticks */
1066 + static AMBVAL   **avlist2;              /* memory positions for sorting */
1067 + static int      i_avlist;               /* index for lists */
1068 +
1069 + static int alatcmp(const void *av1, const void *av2);
1070 +
1071 + static void
1072 + avfree(AMBVAL *av)
1073 + {
1074 +        free(av);
1075 + }
1076 +
1077 + static void
1078 + av2list(
1079 +        AMBVAL *av
1080 + )
1081 + {
1082 + #ifdef DEBUG
1083 +        if (i_avlist >= nambvals)
1084 +                error(CONSISTENCY, "too many ambient values in av2list1");
1085 + #endif
1086 +        avlist1[i_avlist].p = avlist2[i_avlist] = (AMBVAL*)av;
1087 +        avlist1[i_avlist++].t = av->latick;
1088 + }
1089 +
1090 +
1091 + static int
1092 + alatcmp(                        /* compare ambient values for MRA */
1093 +        const void *av1,
1094 +        const void *av2
1095 + )
1096 + {
1097 +        long  lc = ((struct avl *)av2)->t - ((struct avl *)av1)->t;
1098 +        return(lc<0 ? -1 : lc>0 ? 1 : 0);
1099 + }
1100 +
1101 +
1102 + /* GW NOTE 2002/10/3:
1103 + * I used to compare AMBVAL pointers, but found that this was the
1104 + * cause of a serious consistency error with gcc, since the optimizer
1105 + * uses some dangerous trick in pointer subtraction that
1106 + * assumes pointers differ by exact struct size increments.
1107 + */
1108 + static int
1109 + aposcmp(                        /* compare ambient value positions */
1110 +        const void      *avp1,
1111 +        const void      *avp2
1112 + )
1113 + {
1114 +        long    diff = *(char * const *)avp1 - *(char * const *)avp2;
1115 +        if (diff < 0)
1116 +                return(-1);
1117 +        return(diff > 0);
1118 + }
1119 +
1120 +
1121 + static int
1122 + avlmemi(                                /* find list position from address */
1123 +        AMBVAL  *avaddr
1124 + )
1125 + {
1126 +        AMBVAL  **avlpp;
1127 +
1128 +        avlpp = (AMBVAL **)bsearch((char *)&avaddr, (char *)avlist2,
1129 +                        nambvals, sizeof(AMBVAL *), &aposcmp);
1130 +        if (avlpp == NULL)
1131 +                error(CONSISTENCY, "address not found in avlmemi");
1132 +        return(avlpp - avlist2);
1133 + }
1134 +
1135 +
1136 + static void
1137 + sortambvals(                    /* resort ambient values */
1138 +        int     always
1139 + )
1140 + {
1141 +        AMBTREE  oldatrunk;
1142 +        AMBVAL  tav, *tap, *pnext;
1143 +        int     i, j;
1144 +                                        /* see if it's time yet */
1145 +        if (!always && (ambclock++ < lastsort+sortintvl ||
1146 +                        nambvals < SORT_THRESH))
1147 +                return;
1148 +        /*
1149 +         * The idea here is to minimize memory thrashing
1150 +         * in VM systems by improving reference locality.
1151 +         * We do this by periodically sorting our stored ambient
1152 +         * values in memory in order of most recently to least
1153 +         * recently accessed.  This ordering was chosen so that new
1154 +         * ambient values (which tend to be less important) go into
1155 +         * higher memory with the infrequently accessed values.
1156 +         *      Since we expect our values to need sorting less
1157 +         * frequently as the process continues, we double our
1158 +         * waiting interval after each call.
1159 +         *      This routine is also called by setambacc() with
1160 +         * the "always" parameter set to 1 so that the ambient
1161 +         * tree will be rebuilt with the new accuracy parameter.
1162 +         */
1163 +        if (tracktime) {                /* allocate pointer arrays to sort */
1164 +                avlist2 = (AMBVAL **)malloc(nambvals*sizeof(AMBVAL *));
1165 +                avlist1 = (struct avl *)malloc(nambvals*sizeof(struct avl));
1166 +        } else {
1167 +                avlist2 = NULL;
1168 +                avlist1 = NULL;
1169 +        }
1170 +        if (avlist1 == NULL) {          /* no time tracking -- rebuild tree? */
1171 +                if (avlist2 != NULL)
1172 +                        free((void *)avlist2);
1173 +                if (always) {           /* rebuild without sorting */
1174 +                        oldatrunk = atrunk;
1175 +                        atrunk.alist = NULL;
1176 +                        atrunk.kid = NULL;
1177 +                        unloadatree(&oldatrunk, &avinsert);
1178 +                }
1179 +        } else {                        /* sort memory by last access time */
1180 +                /*
1181 +                 * Sorting memory is tricky because it isn't contiguous.
1182 +                 * We have to sort an array of pointers by MRA and also
1183 +                 * by memory position.  We then copy values in "loops"
1184 +                 * to minimize memory hits.  Nevertheless, we will visit
1185 +                 * everyone at least twice, and this is an expensive process
1186 +                 * when we're thrashing, which is when we need to do it.
1187 +                 */
1188 + #ifdef DEBUG
1189 +                sprintf(errmsg, "sorting %u ambient values at ambclock=%lu...",
1190 +                                nambvals, ambclock);
1191 +                eputs(errmsg);
1192 + #endif
1193 +                i_avlist = 0;
1194 +                unloadatree(&atrunk, &av2list); /* empty current tree */
1195 + #ifdef DEBUG
1196 +                if (i_avlist < nambvals)
1197 +                        error(CONSISTENCY, "missing ambient values in sortambvals");
1198 + #endif
1199 +                qsort((char *)avlist1, nambvals, sizeof(struct avl), alatcmp);
1200 +                qsort((char *)avlist2, nambvals, sizeof(AMBVAL *), aposcmp);
1201 +                for (i = 0; i < nambvals; i++) {
1202 +                        if (avlist1[i].p == NULL)
1203 +                                continue;
1204 +                        tap = avlist2[i];
1205 +                        tav = *tap;
1206 +                        for (j = i; (pnext = avlist1[j].p) != tap;
1207 +                                        j = avlmemi(pnext)) {
1208 +                                *(avlist2[j]) = *pnext;
1209 +                                avinsert(avlist2[j]);
1210 +                                avlist1[j].p = NULL;
1211 +                        }
1212 +                        *(avlist2[j]) = tav;
1213 +                        avinsert(avlist2[j]);
1214 +                        avlist1[j].p = NULL;
1215 +                }
1216 +                free((void *)avlist1);
1217 +                free((void *)avlist2);
1218 +                                                /* compute new sort interval */
1219 +                sortintvl = ambclock - lastsort;
1220 +                if (sortintvl >= MAX_SORT_INTVL/2)
1221 +                        sortintvl = MAX_SORT_INTVL;
1222 +                else
1223 +                        sortintvl <<= 1;        /* wait twice as long next */
1224 + #ifdef DEBUG
1225 +                eputs("done\n");
1226 + #endif
1227 +        }
1228 +        if (ambclock >= MAXACLOCK)
1229 +                ambclock = MAXACLOCK/2;
1230 +        lastsort = ambclock;
1231 + }
1232 +
1233 +
1234 + #ifdef  F_SETLKW
1235 +
1236 + static void
1237 + aflock(                 /* lock/unlock ambient file */
1238 +        int  typ
1239 + )
1240 + {
1241 +        static struct flock  fls;       /* static so initialized to zeroes */
1242 +
1243 +        if (typ == fls.l_type)          /* already called? */
1244 +                return;
1245 +        fls.l_type = typ;
1246 +        if (fcntl(fileno(ambfp), F_SETLKW, &fls) < 0)
1247 +                error(SYSTEM, "cannot (un)lock ambient file");
1248 + }
1249 +
1250 +
1251 + int
1252 + ambsync(void)                   /* synchronize ambient file */
1253 + {
1254 +        long  flen;
1255 +        AMBVAL  avs;
1256 +        int  n;
1257 +
1258 +        if (ambfp == NULL)      /* no ambient file? */
1259 +                return(0);
1260 +                                /* gain appropriate access */
1261 +        aflock(nunflshed ? F_WRLCK : F_RDLCK);
1262 +                                /* see if file has grown */
1263 +        if ((flen = lseek(fileno(ambfp), (off_t)0, SEEK_END)) < 0)
1264 +                goto seekerr;
1265 +        if ((n = flen - lastpos) > 0) {         /* file has grown */
1266 +                if (ambinp == NULL) {           /* use duplicate filedes */
1267 +                        ambinp = fdopen(dup(fileno(ambfp)), "r");
1268 +                        if (ambinp == NULL)
1269 +                                error(SYSTEM, "fdopen failed in ambsync");
1270 +                }
1271 +                if (fseek(ambinp, lastpos, SEEK_SET) < 0)
1272 +                        goto seekerr;
1273 +                while (n >= AMBVALSIZ) {        /* load contributed values */
1274 +                        if (!readambval(&avs, ambinp)) {
1275 +                                sprintf(errmsg,
1276 +                        "ambient file \"%s\" corrupted near character %ld",
1277 +                                                ambfile, flen - n);
1278 +                                error(WARNING, errmsg);
1279 +                                break;
1280 +                        }
1281 +                        avstore(&avs);
1282 +                        n -= AMBVALSIZ;
1283 +                }
1284 +                lastpos = flen - n;
1285 +                /*** seek always as safety measure
1286 +                if (n) ***/                     /* alignment */
1287 +                        if (lseek(fileno(ambfp), (off_t)lastpos, SEEK_SET) < 0)
1288 +                                goto seekerr;
1289 +        }
1290 +        n = fflush(ambfp);                      /* calls write() at last */
1291 +        if (n != EOF)
1292 +                lastpos += (long)nunflshed*AMBVALSIZ;
1293 +        else if ((lastpos = lseek(fileno(ambfp), (off_t)0, SEEK_CUR)) < 0)
1294 +                goto seekerr;
1295 +                
1296 +        aflock(F_UNLCK);                        /* release file */
1297 +        nunflshed = 0;
1298 +        return(n);
1299 + seekerr:
1300 +        error(SYSTEM, "seek failed in ambsync");
1301 +        return -1; /* pro forma return */
1302 + }
1303 +
1304 + #else   /* ! F_SETLKW */
1305 +
1306 + int
1307 + ambsync(void)                   /* flush ambient file */
1308 + {
1309 +        if (ambfp == NULL)
1310 +                return(0);
1311 +        nunflshed = 0;
1312 +        return(fflush(ambfp));
1313 + }
1314 +
1315 + #endif  /* ! F_SETLKW */

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines