ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/rt/ambient.c
Revision: 2.131
Committed: Thu Jan 23 22:03:08 2025 UTC (3 months, 1 week ago) by greg
Content type: text/plain
Branch: MAIN
CVS Tags: HEAD
Changes since 2.130: +109 -75 lines
Log Message:
fix: Reverted to tedious file locking, since faster method proved unreliable

File Contents

# Content
1 static const char RCSid[] = "$Id: ambient.c,v 2.128 2024/12/19 19:34:31 greg Exp $";
2 /*
3 * ambient.c - routines dealing with ambient (inter-reflected) component.
4 *
5 * Declarations of external symbols in ambient.h
6 */
7
8 #include "copyright.h"
9
10 #include <string.h>
11
12 #include "platform.h"
13 #include "ray.h"
14 #include "otypes.h"
15 #include "otspecial.h"
16 #include "resolu.h"
17 #include "ambient.h"
18 #include "random.h"
19 #include "pmapamb.h"
20
21 #ifndef OCTSCALE
22 #define OCTSCALE 1.0 /* ceil((valid rad.)/(cube size)) */
23 #endif
24
25 #ifndef MAXASET
26 #define MAXASET 4095 /* maximum number of elements in ambient set */
27 #endif
28 OBJECT ambset[MAXASET+1]={0}; /* ambient include/exclude set */
29
30 double maxarad; /* maximum ambient radius */
31 double minarad; /* minimum ambient radius */
32
33 static AMBTREE atrunk; /* our ambient trunk node */
34
35 static FILE *ambfp = NULL; /* ambient file pointer */
36 static int nunflshed; /* number of unflushed ambient values */
37
38 static double avsum = 0.; /* computed ambient value sum (log) */
39 static unsigned int navsum = 0; /* number of values in avsum */
40 static unsigned int nambvals = 0; /* total number of indirect values */
41 static FILE *ambinp = NULL; /* auxiliary file for input */
42 static off_t lastpos = -1; /* last flush position */
43
44 #define AMBFLUSH (BUFSIZ/AMBVALSIZ)
45
46 #define AVSIZE (sizeof(AMBVAL)-sizeof(SCOLOR)+sizeof(COLORV)*NCSAMP)
47 #define newambval() (AMBVAL *)malloc(AVSIZE)
48
49 #define tfunc(x0, x, x1) (((x)-(x0))/((x1)-(x0)))
50
51 static void initambfile(int cre8);
52 static void avsave(AMBVAL *av);
53 static AMBVAL *avstore(AMBVAL *aval);
54 static AMBTREE *newambtree(void);
55 static void freeambtree(AMBTREE *atp);
56
57 typedef void unloadtf_t(AMBVAL *);
58 static unloadtf_t avinsert;
59 static unloadtf_t avfree;
60 static void unloadatree(AMBTREE *at, unloadtf_t *f);
61
62 static void sortambvals(void);
63
64 static int plugaleak(RAY *r, AMBVAL *ap, FVECT anorm, double ang);
65 static double sumambient(SCOLOR acol, RAY *r, FVECT rn, int al,
66 AMBTREE *at, FVECT c0, double s);
67 static int makeambient(SCOLOR acol, RAY *r, FVECT rn, int al);
68 static int extambient(SCOLOR cr, AMBVAL *ap, FVECT pv, FVECT nv,
69 FVECT uvw[3]);
70
71 #ifdef F_SETLKW
72 static void aflock(int typ);
73 #endif
74
75
76 void
77 setambres( /* set ambient resolution */
78 int ar
79 )
80 {
81 ambres = ar < 0 ? 0 : ar; /* may be done already */
82 /* set min & max radii */
83 if (ar <= 0) {
84 minarad = 0;
85 maxarad = thescene.cusize*0.2;
86 } else {
87 minarad = thescene.cusize / ar;
88 maxarad = 64.0 * minarad; /* heuristic */
89 if (maxarad > thescene.cusize*0.2)
90 maxarad = thescene.cusize*0.2;
91 }
92 if (minarad <= FTINY)
93 minarad = 10.0*FTINY;
94 if (maxarad <= minarad)
95 maxarad = 64.0 * minarad;
96 }
97
98
99 void
100 setambacc( /* set ambient accuracy */
101 double newa
102 )
103 {
104 static double olda; /* remember previous setting here */
105
106 newa *= (newa > 0);
107 if (fabs(newa - olda) >= .05*(newa + olda)) {
108 ambacc = newa;
109 if (ambacc > FTINY && nambvals > 0)
110 sortambvals(); /* rebuild tree */
111 }
112 }
113
114
115 void
116 setambient(void) /* initialize calculation */
117 {
118 int readonly = 0;
119 off_t flen;
120 AMBVAL amb;
121 /* make sure we're fresh */
122 ambdone();
123 /* init ambient limits */
124 setambres(ambres);
125 setambacc(ambacc);
126 if (ambfile == NULL || !ambfile[0])
127 return;
128 if (ambacc <= FTINY) {
129 sprintf(errmsg, "zero ambient accuracy so \"%s\" not opened",
130 ambfile);
131 error(WARNING, errmsg);
132 return;
133 }
134 /* open ambient file */
135 if ((ambfp = fopen(ambfile, "r+")) == NULL)
136 readonly = (ambfp = fopen(ambfile, "r")) != NULL;
137 if (ambfp != NULL) {
138 initambfile(0); /* file exists */
139 lastpos = ftello(ambfp);
140 while (readambval(&amb, ambfp))
141 avstore(&amb);
142 if (readonly) {
143 sprintf(errmsg,
144 "loaded %u values from read-only ambient file",
145 nambvals);
146 error(WARNING, errmsg);
147 fclose(ambfp); /* close file so no writes */
148 ambfp = NULL;
149 return; /* avoid ambsync() */
150 }
151 /* align file pointer */
152 lastpos += (off_t)nambvals*AMBVALSIZ;
153 flen = lseek(fileno(ambfp), 0, SEEK_END);
154 if (flen != lastpos) {
155 sprintf(errmsg,
156 "ignoring last %ld values in ambient file (corrupted)",
157 (flen - lastpos)/AMBVALSIZ);
158 error(WARNING, errmsg);
159 fseeko(ambfp, lastpos, SEEK_SET);
160 ftruncate(fileno(ambfp), lastpos);
161 }
162 } else if ((ambfp = fopen(ambfile, "w+")) != NULL) {
163 initambfile(1); /* else create new file */
164 fflush(ambfp);
165 lastpos = ftello(ambfp);
166 } else {
167 sprintf(errmsg, "cannot open ambient file \"%s\"", ambfile);
168 error(SYSTEM, errmsg);
169 }
170 #ifdef F_SETLKW
171 aflock(F_UNLCK); /* release file */
172 #endif
173 }
174
175
176 void
177 ambdone(void) /* close ambient file and free memory */
178 {
179 if (ambfp != NULL) { /* close ambient file */
180 ambsync();
181 fclose(ambfp);
182 ambfp = NULL;
183 if (ambinp != NULL) {
184 fclose(ambinp);
185 ambinp = NULL;
186 }
187 lastpos = -1;
188 }
189 /* free ambient tree */
190 unloadatree(&atrunk, avfree);
191 /* reset state variables */
192 avsum = 0.;
193 navsum = 0;
194 nambvals = 0;
195 }
196
197
198 void
199 ambnotify( /* record new modifier */
200 OBJECT obj
201 )
202 {
203 static int hitlimit = 0;
204 OBJREC *o;
205 char **amblp;
206
207 if (obj == OVOID) { /* starting over */
208 ambset[0] = 0;
209 hitlimit = 0;
210 return;
211 }
212 o = objptr(obj);
213 if (hitlimit || !ismodifier(o->otype))
214 return;
215 for (amblp = amblist; *amblp != NULL; amblp++)
216 if (!strcmp(o->oname, *amblp)) {
217 if (ambset[0] >= MAXASET) {
218 error(WARNING, "too many modifiers in ambient list");
219 hitlimit++;
220 return; /* should this be fatal? */
221 }
222 insertelem(ambset, obj);
223 return;
224 }
225 }
226
227
228 void
229 multambient( /* compute ambient component & multiply by coef. */
230 SCOLOR aval,
231 RAY *r,
232 FVECT nrm
233 )
234 {
235 static double logAvgAbsorp = 1;
236 static int rdepth = 0; /* ambient recursion */
237 SCOLOR acol, caustic;
238 int i, ok;
239 double d, l;
240
241 /* PMAP: Factor in ambient from photon map, if enabled and ray is
242 * ambient. Return as all ambient components accounted for, else
243 * continue. */
244 if (ambPmap(aval, r, rdepth))
245 return;
246
247 if (logAvgAbsorp > 0) /* exclude in -aw to avoid growth */
248 logAvgAbsorp = log(1.-AVGREFL);
249
250 /* PMAP: Factor in specular-diffuse ambient (caustics) from photon
251 * map, if enabled and ray is primary, else caustic is zero. Continue
252 * with RADIANCE ambient calculation */
253 {/* XXX TEMPORARY */
254 COLOR pmc;
255 scolor_color(pmc, aval);
256 ambPmapCaustic(pmc, r, rdepth);
257 setscolor(caustic, colval(pmc,RED), colval(pmc,GRN), colval(pmc,BLU));
258 }
259 if (ambdiv <= 0) /* no ambient calculation */
260 goto dumbamb;
261 /* check number of bounces */
262 if (rdepth >= ambounce)
263 goto dumbamb;
264 /* check ambient list */
265 if (ambincl != -1 && r->ro != NULL &&
266 ambincl != inset(ambset, r->ro->omod))
267 goto dumbamb;
268
269 if (ambacc <= FTINY) { /* no ambient storage? */
270 double rdot = DOT(nrm,r->ron);
271 int sgn = 1 - 2*(rdot < 0);
272 float dgrad[2], *dgp = NULL;
273 FVECT uvd[2];
274
275 if (sgn*rdot < 0.9999)
276 dgp = dgrad; /* compute rotational grad. */
277 copyscolor(acol, aval);
278 rdepth++;
279 ok = doambient(acol, r, r->rweight*sgn,
280 uvd, NULL, NULL, dgp, NULL);
281 rdepth--;
282 if (!ok)
283 goto dumbamb;
284 if ((ok > 0) & (dgp != NULL)) { /* apply texture */
285 FVECT v1;
286 VCROSS(v1, r->ron, nrm);
287 d = 1.0;
288 for (i = 3; i--; )
289 d += sgn*v1[i] * (dgp[0]*uvd[0][i] + dgp[1]*uvd[1][i]);
290 if (d >= 0.05)
291 scalescolor(acol, d);
292 }
293 copyscolor(aval, acol);
294
295 /* PMAP: add in caustic */
296 saddscolor(aval, caustic);
297 return;
298 }
299 /* interpolate ambient value */
300 scolorblack(acol);
301 d = sumambient(acol, r, nrm, rdepth,
302 &atrunk, thescene.cuorg, thescene.cusize);
303
304 if (d > FTINY) {
305 scalescolor(acol, 1.0/d);
306 smultscolor(aval, acol);
307
308 /* PMAP: add in caustic */
309 saddscolor(aval, caustic);
310 return;
311 }
312
313 rdepth++; /* need to cache new value */
314 ok = makeambient(acol, r, nrm, rdepth-1);
315 rdepth--;
316
317 if (ok) {
318 smultscolor(aval, acol); /* computed new value */
319
320 /* PMAP: add in caustic */
321 saddscolor(aval, caustic);
322 return;
323 }
324
325 dumbamb: /* return global value */
326 if ((ambvwt <= 0) | (navsum == 0)) {
327 smultcolor(aval, ambval);
328
329 /* PMAP: add in caustic */
330 saddscolor(aval, caustic);
331 return;
332 }
333
334 l = bright(ambval); /* average in computations */
335 if (l > FTINY) {
336 d = (log(l)*(double)ambvwt + avsum + logAvgAbsorp*navsum) /
337 (double)(ambvwt + navsum);
338 d = exp(d) / l;
339 scalescolor(aval, d);
340 smultcolor(aval, ambval); /* apply color of ambval */
341 } else {
342 d = exp( avsum/(double)navsum + logAvgAbsorp );
343 scalescolor(aval, d); /* neutral color */
344 }
345 }
346
347
348 /* Plug a potential leak where ambient cache value is occluded */
349 static int
350 plugaleak(RAY *r, AMBVAL *ap, FVECT anorm, double ang)
351 {
352 const double cost70sq = 0.1169778; /* cos(70deg)^2 */
353 RAY rtst;
354 FVECT vdif;
355 double normdot, ndotd, nadotd;
356 double a, b, c, t[2];
357
358 ang += 2.*PI*(ang < 0); /* check direction flags */
359 if ( !(ap->corral>>(int)(ang*(16./PI)) & 1) )
360 return(0);
361 /*
362 * Generate test ray, targeting 20 degrees above sample point plane
363 * along surface normal from cache position. This should be high
364 * enough to miss local geometry we don't really care about.
365 */
366 VSUB(vdif, ap->pos, r->rop);
367 normdot = DOT(anorm, r->ron);
368 ndotd = DOT(vdif, r->ron);
369 nadotd = DOT(vdif, anorm);
370 a = normdot*normdot - cost70sq;
371 b = 2.0*(normdot*ndotd - nadotd*cost70sq);
372 c = ndotd*ndotd - DOT(vdif,vdif)*cost70sq;
373 if (quadratic(t, a, b, c) != 2)
374 return(1); /* should rarely happen */
375 if (t[1] <= FTINY)
376 return(0); /* should fail behind test */
377 rayorigin(&rtst, SHADOW, r, NULL);
378 VSUM(rtst.rdir, vdif, anorm, t[1]); /* further dist. > plane */
379 rtst.rmax = normalize(rtst.rdir); /* short ray test */
380 while (localhit(&rtst, &thescene)) { /* check for occluder */
381 OBJREC *m = findmaterial(rtst.ro);
382 if (m != NULL && !istransp(m) && !isBSDFproxy(m) &&
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 SCOLOR 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 const int sgn = 1 - 2*(DOT(r->ron,rn) < 0);
407 double maxangle = minangle + ambacc;
408 double wsum = 0.0;
409 FVECT ck0;
410 int i, j;
411 AMBVAL *av;
412
413 if (at->kid != NULL) { /* sum children first */
414 s *= 0.5;
415 for (i = 0; i < 8; i++) {
416 for (j = 0; j < 3; j++) {
417 ck0[j] = c0[j];
418 if (1<<j & i)
419 ck0[j] += s;
420 if (r->rop[j] < ck0[j] - OCTSCALE*s)
421 break;
422 if (r->rop[j] > ck0[j] + (1.0+OCTSCALE)*s)
423 break;
424 }
425 if (j == 3)
426 wsum += sumambient(acol, r, rn, al,
427 at->kid+i, ck0, s);
428 }
429 /* good enough? */
430 if ((wsum >= 0.05) & (s*ambacc > minarad))
431 return(wsum);
432 }
433 /* adjust maximum angle */
434 if (at->alist != NULL && (at->alist->lvl <= al) & (r->rweight < 0.6))
435 maxangle = (maxangle - PI/2.)*pow(r->rweight,0.13) + PI/2.;
436 /* sum this node */
437 for (av = at->alist; av != NULL; av = av->next) {
438 double u, v, d, delta_r2, delta_t2;
439 SCOLOR sct;
440 FVECT uvw[3];
441 /*
442 * Ambient level test
443 */
444 if (av->lvl > al || /* list sorted, so this works */
445 (av->lvl == al) & (av->weight < 0.9*r->rweight))
446 break;
447 /*
448 * Direction test using unperturbed normal
449 */
450 decodedir(uvw[2], av->ndir);
451 d = sgn * DOT(uvw[2], r->ron);
452 if (d <= 0.0) /* >= 90 degrees */
453 continue;
454 delta_r2 = 2.0 - 2.0*d; /* approx. radians^2 */
455 if (delta_r2 >= maxangle*maxangle)
456 continue;
457 /*
458 * Modified ray behind test
459 */
460 VSUB(ck0, r->rop, av->pos);
461 d = DOT(ck0, uvw[2]);
462 if (d < -minarad*ambacc)
463 continue;
464 d /= av->rad[0];
465 delta_t2 = d*d;
466 if (delta_t2 >= ambacc*ambacc)
467 continue;
468 /*
469 * Elliptical radii test based on Hessian
470 */
471 decodedir(uvw[0], av->udir);
472 VCROSS(uvw[1], uvw[2], uvw[0]);
473 d = (u = DOT(ck0, uvw[0])) / av->rad[0];
474 delta_t2 += d*d;
475 d = (v = DOT(ck0, uvw[1])) / av->rad[1];
476 delta_t2 += d*d;
477 if (delta_t2 >= ambacc*ambacc)
478 continue;
479 /*
480 * Test for potential light leak
481 */
482 if (av->corral && plugaleak(r, av, uvw[2], atan2a(v,u)))
483 continue;
484 /*
485 * Extrapolate value and compute final weight (hat function)
486 */
487 if (!extambient(sct, av, r->rop, rn, uvw))
488 continue;
489 d = tfunc(maxangle, sqrt(delta_r2), 0.0) *
490 tfunc(ambacc, sqrt(delta_t2), 0.0);
491 scalescolor(sct, d);
492 saddscolor(acol, sct);
493 wsum += d;
494 }
495 return(wsum);
496 }
497
498
499 static int
500 makeambient( /* make a new ambient value for storage */
501 SCOLOR acol,
502 RAY *r,
503 FVECT rn,
504 int al
505 )
506 {
507 int sgn = 1 - 2*(DOT(r->ron,rn) < 0);
508 AMBVAL amb;
509 FVECT uvw[3];
510 int i;
511
512 amb.weight = 1.0; /* compute weight */
513 for (i = al; i-- > 0; )
514 amb.weight *= AVGREFL;
515 if (r->rweight < 0.1*amb.weight) /* heuristic override */
516 amb.weight = 1.25*r->rweight;
517 setscolor(acol, AVGREFL, AVGREFL, AVGREFL);
518 /* compute ambient */
519 i = doambient(acol, r, amb.weight*sgn,
520 uvw, amb.rad, amb.gpos, amb.gdir, &amb.corral);
521 scalescolor(acol, 1./AVGREFL); /* undo assumed reflectance */
522 if (i <= 0 || amb.rad[0] <= FTINY) /* no Hessian or zero radius */
523 return(i);
524 uvw[2][0] = sgn*r->ron[0]; /* orient unperturbed normal */
525 uvw[2][1] = sgn*r->ron[1];
526 uvw[2][2] = sgn*r->ron[2];
527 /* store value */
528 VCOPY(amb.pos, r->rop);
529 amb.ndir = encodedir(uvw[2]);
530 amb.udir = encodedir(uvw[0]);
531 amb.lvl = al;
532 copyscolor(amb.val, acol);
533 avsave(&amb); /* insert and save to file */
534 if (DOT(uvw[2],rn) < 0.9999) /* texture? */
535 extambient(acol, &amb, r->rop, rn, uvw);
536 return(1);
537 }
538
539
540 static int
541 extambient( /* extrapolate value at pv, nv */
542 SCOLOR scr,
543 AMBVAL *ap,
544 FVECT pv,
545 FVECT nv,
546 FVECT uvw[3]
547 )
548 {
549 const double min_d = 0.05;
550 const double max_d = 20.;
551 static FVECT my_uvw[3];
552 FVECT v1;
553 int i;
554 double d = 1.0; /* zeroeth order */
555
556 if (uvw == NULL) { /* need local coordinates? */
557 decodedir(my_uvw[2], ap->ndir);
558 decodedir(my_uvw[0], ap->udir);
559 VCROSS(my_uvw[1], my_uvw[2], my_uvw[0]);
560 uvw = my_uvw;
561 }
562 for (i = 3; i--; ) /* gradient due to translation */
563 d += (pv[i] - ap->pos[i]) *
564 (ap->gpos[0]*uvw[0][i] + ap->gpos[1]*uvw[1][i]);
565
566 VCROSS(v1, uvw[2], nv); /* gradient due to rotation */
567 for (i = 3; i--; )
568 d += v1[i] * (ap->gdir[0]*uvw[0][i] + ap->gdir[1]*uvw[1][i]);
569
570 if (d < min_d) /* clamp min/max scaling */
571 d = min_d;
572 else if (d > max_d)
573 d = max_d;
574 copyscolor(scr, ap->val);
575 scalescolor(scr, d);
576 return(d > min_d);
577 }
578
579
580 static void
581 avinsert( /* insert ambient value in our tree */
582 AMBVAL *av
583 )
584 {
585 AMBTREE *at;
586 AMBVAL *ap;
587 AMBVAL avh;
588 FVECT ck0;
589 double s;
590 int branch;
591 int i;
592
593 if (av->rad[0] <= FTINY)
594 error(CONSISTENCY, "zero ambient radius in avinsert");
595 at = &atrunk;
596 VCOPY(ck0, thescene.cuorg);
597 s = thescene.cusize;
598 while (s*(OCTSCALE/2) > av->rad[1]*ambacc) {
599 if (at->kid == NULL)
600 if ((at->kid = newambtree()) == NULL)
601 error(SYSTEM, "out of memory in avinsert");
602 s *= 0.5;
603 branch = 0;
604 for (i = 0; i < 3; i++)
605 if (av->pos[i] > ck0[i] + s) {
606 ck0[i] += s;
607 branch |= 1 << i;
608 }
609 at = at->kid + branch;
610 }
611 avh.next = at->alist; /* order by increasing level */
612 for (ap = &avh; ap->next != NULL; ap = ap->next)
613 if ( ap->next->lvl > av->lvl ||
614 (ap->next->lvl == av->lvl) &
615 (ap->next->weight <= av->weight) )
616 break;
617 av->next = ap->next;
618 ap->next = (AMBVAL*)av;
619 at->alist = avh.next;
620 }
621
622
623 static void
624 initambfile( /* initialize ambient file */
625 int cre8
626 )
627 {
628 extern char *progname, *octname;
629 static char *mybuf = NULL;
630 int ntries = 3;
631
632 if (!AMBFLUSH)
633 error(INTERNAL, "BUFSIZ too small in initambfile");
634 #ifdef F_SETLKW
635 aflock(cre8 ? F_WRLCK : F_RDLCK);
636 #endif
637 SET_FILE_BINARY(ambfp);
638 if (mybuf == NULL)
639 mybuf = (char *)bmalloc(BUFSIZ+8);
640 setbuf(ambfp, mybuf);
641 nunflshed = 0;
642 retry:
643 if (cre8) { /* new file */
644 newheader("RADIANCE", ambfp);
645 fprintf(ambfp, "%s -av %g %g %g -aw %d -ab %d -aa %g ",
646 progname, colval(ambval,RED),
647 colval(ambval,GRN), colval(ambval,BLU),
648 ambvwt, ambounce, ambacc);
649 fprintf(ambfp, "-ad %d -as %d -ar %d ",
650 ambdiv, ambssamp, ambres);
651 fprintf(ambfp, "-dr %d -ds %g -dt %g -dc %g ", directrelay,
652 srcsizerat, shadthresh, shadcert);
653 fprintf(ambfp, "-ss %g -st %g -lr %d -lw %g ", specjitter,
654 specthresh, maxdepth, minweight);
655 fprintf(ambfp, "-cw %g %g -cs %d ", WLPART[3], WLPART[0], NCSAMP);
656 if (octname != NULL)
657 fputs(octname, ambfp);
658 fputc('\n', ambfp); /* end of command line, not header! */
659 fprintf(ambfp, "SOFTWARE= %s\n", VersionID);
660 fputnow(ambfp);
661 AMB_CNDX = CNDX; /* use current spectral sampling */
662 AMB_WLPART = WLPART;
663 fputwlsplit(WLPART, ambfp);
664 fputncomp(NCSAMP, ambfp);
665 fputformat(AMBFMT, ambfp);
666 fputc('\n', ambfp);
667 putambmagic(ambfp);
668 } else if (getheader(ambfp, amb_headline, NULL) < 0 || !hasambmagic(ambfp)) {
669 if (--ntries > 0 && ftell(ambfp) == 0) {
670 #ifdef F_SETLKW
671 aflock(F_UNLCK);
672 clearerr(ambfp);
673 sleep(2);
674 aflock(F_RDLCK);
675 #else
676 clearerr(ambfp);
677 sleep(2);
678 #endif
679 goto retry;
680 }
681 error(USER, "bad/incompatible ambient file");
682 }
683 if ((AMB_CNDX != CNDX) | (AMB_WLPART != WLPART)) {
684 if (setspectrsamp(AMB_CNDX, AMB_WLPART) < 0)
685 error(USER, "bad wavelength sampling in ambient file");
686 if (AMB_CNDX[3] == CNDX[3] && FABSEQ(AMB_WLPART[0],WLPART[0]) &&
687 FABSEQ(AMB_WLPART[3],WLPART[3])) {
688 AMB_CNDX = CNDX;
689 AMB_WLPART = WLPART; /* just the same */
690 } else
691 error(WARNING, "different ambient file wavelength sampling");
692 }
693 }
694
695
696 static void
697 avsave( /* insert and save an ambient value */
698 AMBVAL *av
699 )
700 {
701 avstore(av);
702 if (ambfp == NULL)
703 return;
704 if (writambval(av, ambfp) < 0)
705 goto writerr;
706 if (++nunflshed >= AMBFLUSH)
707 if (ambsync() == EOF)
708 goto writerr;
709 return;
710 writerr:
711 error(SYSTEM, "error writing to ambient file");
712 }
713
714
715 static AMBVAL *
716 avstore( /* allocate memory and save aval */
717 AMBVAL *aval
718 )
719 {
720 AMBVAL *av;
721 double d;
722
723 if ((av = newambval()) == NULL)
724 error(SYSTEM, "out of memory in avstore");
725 memcpy(av, aval, AVSIZE); /* AVSIZE <= sizeof(AMBVAL) */
726 av->next = NULL;
727 nambvals++;
728 d = pbright(av->val);
729 if (d > FTINY) { /* add to log sum for averaging */
730 avsum += log(d);
731 navsum++;
732 }
733 avinsert(av); /* insert in our cache tree */
734 return(av);
735 }
736
737
738 #define ATALLOCSZ 512 /* #/8 trees to allocate at once */
739
740 static AMBTREE *atfreelist = NULL; /* free ambient tree structures */
741
742
743 static AMBTREE *
744 newambtree(void) /* allocate 8 ambient tree structs */
745 {
746 AMBTREE *atp, *upperlim;
747
748 if (atfreelist == NULL) { /* get more nodes */
749 atfreelist = (AMBTREE *)malloc(ATALLOCSZ*8*sizeof(AMBTREE));
750 if (atfreelist == NULL)
751 return(NULL);
752 /* link new free list */
753 upperlim = atfreelist + 8*(ATALLOCSZ-1);
754 for (atp = atfreelist; atp < upperlim; atp += 8)
755 atp->kid = atp + 8;
756 atp->kid = NULL;
757 }
758 atp = atfreelist;
759 atfreelist = atp->kid;
760 memset(atp, 0, 8*sizeof(AMBTREE));
761 return(atp);
762 }
763
764
765 static void
766 freeambtree( /* free 8 ambient tree structs */
767 AMBTREE *atp
768 )
769 {
770 atp->kid = atfreelist;
771 atfreelist = atp;
772 }
773
774
775 static void
776 unloadatree( /* unload an ambient value tree */
777 AMBTREE *at,
778 unloadtf_t *f
779 )
780 {
781 AMBVAL *av;
782 int i;
783 /* transfer values at this node */
784 for (av = at->alist; av != NULL; av = at->alist) {
785 at->alist = av->next;
786 av->next = NULL;
787 (*f)(av);
788 }
789 if (at->kid == NULL)
790 return;
791 for (i = 0; i < 8; i++) /* transfer and free children */
792 unloadatree(at->kid+i, f);
793 freeambtree(at->kid);
794 at->kid = NULL;
795 }
796
797
798 static void
799 avfree(AMBVAL *av)
800 {
801 free(av);
802 }
803
804
805 static void
806 sortambvals(void) /* resort ambient values */
807 {
808 AMBTREE oldatrunk = atrunk;
809
810 atrunk.alist = NULL;
811 atrunk.kid = NULL;
812 unloadatree(&oldatrunk, avinsert);
813 }
814
815
816 #ifdef F_SETLKW
817
818 static void
819 aflock( /* lock/unlock ambient file */
820 int typ
821 )
822 {
823 static struct flock fls; /* static so initialized to zeroes */
824
825 if (typ == fls.l_type) /* already called? */
826 return;
827
828 fls.l_type = typ;
829 do
830 if (fcntl(fileno(ambfp), F_SETLKW, &fls) != -1)
831 return;
832 while (errno == EINTR);
833
834 error(SYSTEM, "cannot (un)lock ambient file");
835 }
836
837
838 int
839 ambsync(void) /* synchronize ambient file */
840 {
841 off_t flen;
842 AMBVAL avs;
843 int n;
844
845 if (ambfp == NULL) /* no ambient file? */
846 return(0);
847 /* gain appropriate access */
848 aflock(nunflshed ? F_WRLCK : F_RDLCK);
849 /* see if file has grown */
850 if ((flen = lseek(fileno(ambfp), 0, SEEK_END)) < 0)
851 goto seekerr;
852 if ((n = flen - lastpos) > 0) { /* file has grown */
853 if (ambinp == NULL) { /* get new file pointer */
854 ambinp = fopen(ambfile, "rb");
855 if (ambinp == NULL)
856 error(SYSTEM, "fopen failed in ambsync");
857 }
858 if (fseeko(ambinp, lastpos, SEEK_SET) < 0)
859 goto seekerr;
860 while (n >= AMBVALSIZ) { /* load contributed values */
861 if (!readambval(&avs, ambinp)) {
862 sprintf(errmsg,
863 "ambient file \"%s\" corrupted near character %ld",
864 ambfile, flen - n);
865 error(WARNING, errmsg);
866 break;
867 }
868 avstore(&avs);
869 n -= AMBVALSIZ;
870 }
871 lastpos = flen - n; /* check alignment */
872 if (n && lseek(fileno(ambfp), lastpos, SEEK_SET) < 0)
873 goto seekerr;
874 }
875 n = fflush(ambfp); /* calls write() at last */
876 lastpos += (off_t)nunflshed*AMBVALSIZ;
877 aflock(F_UNLCK); /* release file */
878 nunflshed = 0;
879 return(n);
880 seekerr:
881 error(SYSTEM, "seek failed in ambsync");
882 return(EOF); /* pro forma return */
883 }
884
885 #else /* ! F_SETLKW */
886
887 int
888 ambsync(void) /* flush ambient file */
889 {
890 if (ambfp == NULL)
891 return(0);
892 nunflshed = 0;
893 return(fflush(ambfp));
894 }
895
896 #endif /* ! F_SETLKW */