ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/util/glaresrc.c
Revision: 2.6
Committed: Mon May 29 16:47:54 2006 UTC (17 years, 10 months ago) by greg
Content type: text/plain
Branch: MAIN
CVS Tags: rad5R2, rad4R2P2, rad5R0, rad5R1, rad4R2, rad4R1, rad4R0, rad3R8, rad3R9, rad4R2P1, rad5R3
Changes since 2.5: +3 -1 lines
Log Message:
Fixed -DDEBUG compile

File Contents

# User Rev Content
1 greg 1.1 #ifndef lint
2 greg 2.6 static const char RCSid[] = "$Id: glaresrc.c,v 2.5 2004/01/02 12:51:54 schorsch Exp $";
3 greg 1.1 #endif
4     /*
5     * Gather samples and compute glare sources.
6     */
7    
8     #include "glare.h"
9 greg 1.15 #include "linregr.h"
10 greg 1.1
11     #define vcont(vd,vu) ((vu)-(vd)<=SEPS)
12     #define hcont(s1,s2) ((s1)->r-(s2)->l>=-SEPS&&(s2)->r-(s1)->l>=-SEPS)
13    
14     struct source *curlist = NULL; /* current source list */
15     struct source *donelist = NULL; /* finished sources */
16    
17 greg 2.6 void pict_stats(void);
18    
19 schorsch 2.5 static struct srcspan * newspan(int l, int r, int v, float *sb);
20     static struct srcspan * newspan(int l, int r, int v, float *sb);
21     static void addindirect(int h, int v, double br);
22     static void addsrcspan(struct srcspan *nss);
23     static void mergesource(struct source *sp, struct source *ap);
24     static void close_sources(int v);
25     static void close_allsrcs(void);
26     static struct srcspan * splitspan(register struct srcspan *sso, double h, double v, double m);
27     static struct source * splitsource(struct source *so);
28     static void donesource(register struct source *sp);
29     static struct source * findbuddy(register struct source *s, register struct source *l);
30     static void absorb(register struct source *s);
31     static void freespans(struct source *sp);
32    
33    
34     static struct srcspan *
35     newspan( /* allocate a new source span */
36     int l,
37     int r,
38     int v,
39     float *sb
40     )
41 greg 1.1 {
42     register struct srcspan *ss;
43     register int i;
44    
45     ss = (struct srcspan *)malloc(sizeof(struct srcspan));
46     if (ss == NULL)
47     memerr("source spans");
48     ss->l = l;
49     ss->r = r;
50     ss->v = v;
51     ss->brsum = 0.0;
52     for (i = l; i < r; i++)
53     ss->brsum += sb[i+hsize];
54     return(ss);
55     }
56    
57    
58 schorsch 2.5 extern void
59     analyze(void) /* analyze our scene */
60 greg 1.1 {
61     int h, v;
62     int left;
63     float *spanbr;
64    
65     spanbr = (float *)malloc((2*hsize+1)*sizeof(float));
66     if (spanbr == NULL)
67     memerr("view span brightness buffer");
68     for (v = vsize; v >= -vsize; v--) {
69 greg 1.5 close_sources(v);
70 greg 1.8 #ifndef DEBUG
71 greg 2.2 if (verbose) {
72 greg 1.8 fprintf(stderr, "%s: analyzing... %3ld%%\r",
73     progname, 100L*(vsize-v)/(2*vsize));
74 greg 2.2 fflush(stderr);
75     }
76 greg 1.8 #endif
77 greg 1.1 getviewspan(v, spanbr);
78     left = hsize + 1;
79 greg 1.2 for (h = -hsize; h <= hsize; h++) {
80 greg 1.1 if (spanbr[h+hsize] < 0.0) { /* off view */
81     if (left < h) {
82 greg 1.2 addsrcspan(newspan(left,h,v,spanbr));
83 greg 1.1 left = hsize + 1;
84     }
85     continue;
86     }
87     if (spanbr[h+hsize] > threshold) { /* in source */
88 greg 1.2 if (left > h)
89 greg 1.1 left = h;
90     } else { /* out of source */
91     if (left < h) {
92 greg 1.2 addsrcspan(newspan(left,h,v,spanbr));
93 greg 1.1 left = hsize + 1;
94     }
95 greg 1.13 addindirect(h, v, spanbr[h+hsize]);
96 greg 1.1 }
97 greg 1.2 }
98 greg 1.1 if (left < h)
99 greg 1.2 addsrcspan(newspan(left,h,v,spanbr));
100 greg 1.1 }
101 greg 2.4 free((void *)spanbr);
102 greg 1.1 close_allsrcs();
103     }
104    
105    
106 schorsch 2.5 static void
107     addindirect( /* add brightness to indirect illuminances */
108     int h,
109     int v,
110     double br
111     )
112 greg 1.1 {
113     double tanb, d;
114 greg 1.13 int hl;
115 greg 1.1 register int i;
116    
117 greg 1.13 hl = hlim(v);
118     if (h <= -hl) { /* left region */
119     d = (double)(-h-hl)/sampdens;
120 greg 1.9 if (d >= 1.0-FTINY)
121 greg 1.1 return;
122     tanb = d/sqrt(1.0-d*d);
123     for (i = 0; i < nglardirs; i++) {
124     d = indirect[i].lcos - tanb*indirect[i].lsin;
125     if (d > 0.0) {
126     indirect[i].sum += d * br;
127 greg 1.10 indirect[i].n += d;
128 greg 1.1 }
129     }
130     return;
131     }
132 greg 1.13 if (h >= hl) { /* right region */
133     d = (double)(-h+hl)/sampdens;
134 greg 1.9 if (d <= -1.0+FTINY)
135 greg 1.1 return;
136     tanb = d/sqrt(1.0-d*d);
137     for (i = 0; i < nglardirs; i++) {
138     d = indirect[i].rcos - tanb*indirect[i].rsin;
139     if (d > 0.0) {
140     indirect[i].sum += d * br;
141 greg 1.10 indirect[i].n += d;
142 greg 1.1 }
143     }
144     return;
145     }
146     /* central region */
147     for (i = 0; i < nglardirs; i++) {
148 greg 1.13 d = cos(h_theta(h,v) - indirect[i].theta);
149 greg 1.1 if (d > 0.0) {
150     indirect[i].sum += d * br;
151 greg 1.10 indirect[i].n += d;
152 greg 1.1 }
153     }
154     }
155    
156    
157 schorsch 2.5 extern void
158     comp_thresh(void) /* compute glare threshold */
159 greg 1.1 {
160     int h, v;
161     int nsamps;
162     double brsum, br;
163    
164     if (verbose)
165     fprintf(stderr, "%s: computing glare threshold...\n",
166     progname);
167     brsum = 0.0;
168     nsamps = 0;
169 greg 1.12 for (v = vsize; v >= -vsize; v -= TSAMPSTEP) {
170 greg 1.1 for (h = -hsize; h <= hsize; h += TSAMPSTEP) {
171 greg 1.13 if ((br = getviewpix(h, v)) < 0.0)
172 greg 1.1 continue;
173     brsum += br;
174     nsamps++;
175     }
176 greg 1.12 }
177 greg 1.1 if (nsamps == 0) {
178     fprintf(stderr, "%s: no viewable scene!\n", progname);
179     exit(1);
180     }
181     threshold = GLAREBR * brsum / nsamps;
182     if (threshold <= FTINY) {
183     fprintf(stderr, "%s: threshold zero!\n", progname);
184     exit(1);
185     }
186 greg 1.2 if (verbose) {
187 greg 1.8 #ifdef DEBUG
188 greg 1.2 pict_stats();
189 greg 1.8 #endif
190 greg 1.1 fprintf(stderr,
191     "%s: threshold set to %f cd/m2 from %d samples\n",
192     progname, threshold, nsamps);
193 greg 1.2 }
194 greg 1.1 }
195    
196    
197 schorsch 2.5 static void
198     addsrcspan( /* add new source span to our list */
199     struct srcspan *nss
200     )
201 greg 1.1 {
202     struct source *last, *cs, *this;
203     register struct srcspan *ss;
204    
205     cs = NULL;
206     for (this = curlist; this != NULL; this = this->next) {
207     for (ss = this->first; ss != NULL; ss = ss->next) {
208     if (!vcont(nss->v, ss->v))
209     break;
210     if (hcont(ss, nss)) {
211     if (cs == NULL)
212     cs = this;
213     else {
214 greg 1.7 last->next = this->next;
215 greg 1.1 mergesource(cs, this);
216 greg 1.5 this = last;
217 greg 1.1 }
218     break;
219     }
220     }
221     last = this;
222     }
223     if (cs == NULL) {
224     cs = (struct source *)malloc(sizeof(struct source));
225     if (cs == NULL)
226     memerr("source records");
227 greg 1.7 cs->dom = 0.0;
228 greg 1.1 cs->first = NULL;
229     cs->next = curlist;
230     curlist = cs;
231     }
232     nss->next = cs->first;
233     cs->first = nss;
234     }
235    
236    
237 schorsch 2.5 static void
238     mergesource( /* merge source ap into source sp */
239     struct source *sp,
240     struct source *ap
241     )
242 greg 1.1 {
243     struct srcspan head;
244     register struct srcspan *alp, *prev, *tp;
245    
246     head.next = sp->first;
247     prev = &head;
248     alp = ap->first;
249     while (alp != NULL && prev->next != NULL) {
250     if (prev->next->v > alp->v) {
251     tp = alp->next;
252     alp->next = prev->next;
253     prev->next = alp;
254     alp = tp;
255     }
256     prev = prev->next;
257     }
258     if (prev->next == NULL)
259     prev->next = alp;
260     sp->first = head.next;
261 greg 1.7 if (ap->dom > 0.0 && sp->dom > 0.0) { /* sources are done */
262     sp->dir[0] *= sp->dom;
263     sp->dir[1] *= sp->dom;
264     sp->dir[2] *= sp->dom;
265     fvsum(sp->dir, sp->dir, ap->dir, ap->dom);
266     normalize(sp->dir);
267     sp->brt = (sp->brt*sp->dom + ap->brt*ap->dom)
268     / (sp->dom + ap->dom);
269     }
270 greg 2.4 free((void *)ap);
271 greg 1.1 }
272    
273    
274 schorsch 2.5 static void
275     close_sources( /* close sources above v */
276     int v
277     )
278 greg 1.1 {
279     struct source head;
280     register struct source *last, *this;
281    
282     head.next = curlist;
283     last = &head;
284     for (this = curlist; this != NULL; this = this->next)
285     if (!vcont(v, this->first->v)) {
286     last->next = this->next;
287     donesource(this);
288     this = last;
289     } else
290     last = this;
291     curlist = head.next;
292     }
293    
294    
295 schorsch 2.5 static void
296     close_allsrcs(void) /* done with everything */
297 greg 1.1 {
298 greg 1.5 register struct source *this, *next;
299 greg 1.1
300 greg 1.5 this = curlist;
301     while (this != NULL) {
302     next = this->next;
303 greg 1.1 donesource(this);
304 greg 1.5 this = next;
305 greg 1.4 }
306 greg 1.1 curlist = NULL;
307     }
308    
309    
310 schorsch 2.5 static struct srcspan *
311     splitspan( /* divide source span at point */
312     register struct srcspan *sso,
313     double h,
314     double v,
315     double m
316     )
317 greg 1.15 {
318     register struct srcspan *ssn;
319     double d;
320     int hs;
321    
322     d = h - m*(sso->v - v);
323     hs = d < 0. ? d-.5 : d+.5;
324     if (sso->l >= hs)
325     return(NULL);
326     if (sso->r <= hs)
327     return(sso);
328     /* need to split it */
329     ssn = (struct srcspan *)malloc(sizeof(struct srcspan));
330     if (ssn == NULL)
331     memerr("source spans in splitspan");
332     ssn->brsum = (double)(hs - sso->l)/(sso->r - sso->l) * sso->brsum;
333     sso->brsum -= ssn->brsum;
334     ssn->v = sso->v;
335     ssn->l = sso->l;
336     ssn->r = sso->l = hs;
337     return(ssn);
338     }
339    
340    
341 schorsch 2.5 static struct source *
342     splitsource( /* divide source in two if it's big and long */
343     struct source *so
344     )
345 greg 1.15 {
346     LRSUM lr;
347     LRLIN fit;
348     register struct srcspan *ss, *ssn;
349     struct srcspan *ssl, *ssnl, head;
350     int h;
351     double mh, mv;
352     struct source *sn;
353    
354     lrclear(&lr);
355     for (ss = so->first; ss != NULL; ss = ss->next)
356     for (h = ss->l; h < ss->r; h++)
357     lrpoint(h, ss->v, &lr);
358     if ((double)lr.n/(sampdens*sampdens) < SABIG)
359     return(NULL); /* too small */
360     if (lrfit(&fit, &lr) < 0)
361     return(NULL); /* can't fit a line */
362     if (fit.correlation < LCORR && fit.correlation > -LCORR)
363     return(NULL);
364     if (verbose)
365     fprintf(stderr, "%s: splitting large source\n", progname);
366     mh = lrxavg(&lr);
367     mv = lryavg(&lr);
368     sn = (struct source *)malloc(sizeof(struct source));
369     if (sn == NULL)
370     memerr("source records in splitsource");
371     sn->dom = 0.0;
372     sn->first = NULL;
373     ssnl = NULL;
374     head.next = so->first;
375     ssl = &head;
376     for (ss = so->first; ss != NULL; ssl = ss, ss = ss->next)
377     if ((ssn = splitspan(ss, mh, mv, fit.slope)) != NULL) {
378     if (ssn == ss) { /* remove from old */
379     ssl->next = ss->next;
380     ss = ssl;
381     }
382     if (ssnl == NULL) /* add to new */
383     sn->first = ssn;
384     else
385     ssnl->next = ssn;
386     ssn->next = NULL;
387     ssnl = ssn;
388     }
389     so->first = head.next;
390     return(sn);
391     }
392    
393    
394 schorsch 2.5 static void
395     donesource( /* finished with this source */
396     register struct source *sp
397     )
398 greg 1.1 {
399 greg 1.15 struct source *newsrc;
400 greg 1.1 register struct srcspan *ss;
401     int h, n;
402 greg 1.13 double hsum, vsum, d;
403 greg 1.1
404 greg 1.15 while ((newsrc = splitsource(sp)) != NULL) /* split it? */
405     donesource(newsrc);
406 greg 1.1 sp->dom = 0.0;
407 greg 1.13 hsum = vsum = 0.0;
408 greg 1.1 sp->brt = 0.0;
409 greg 1.3 n = 0;
410 greg 1.1 for (ss = sp->first; ss != NULL; ss = ss->next) {
411     sp->brt += ss->brsum;
412     n += ss->r - ss->l;
413 greg 1.14 for (h = ss->l; h < ss->r; h++) {
414 greg 1.13 d = pixsize(h, ss->v);
415     hsum += d*h;
416     vsum += d*ss->v;
417     sp->dom += d;
418     }
419 greg 1.1 }
420 greg 2.3 freespans(sp);
421     if (sp->dom <= FTINY) { /* must be right at edge of image */
422 greg 2.4 free((void *)sp);
423 greg 2.3 return;
424     }
425 greg 1.1 sp->brt /= (double)n;
426 greg 1.13 compdir(sp->dir, (int)(hsum/sp->dom), (int)(vsum/sp->dom));
427 greg 1.1 sp->next = donelist;
428     donelist = sp;
429     if (verbose)
430     fprintf(stderr,
431 greg 1.12 "%s: source at [%.3f,%.3f,%.3f], dw %.5f, br %.1f (%d samps)\n",
432 greg 1.2 progname, sp->dir[0], sp->dir[1], sp->dir[2],
433 greg 1.5 sp->dom, sp->brt, n);
434 greg 1.7 }
435    
436    
437 schorsch 2.5 static struct source *
438     findbuddy( /* find close enough source to s in l*/
439     register struct source *s,
440     register struct source *l
441     )
442 greg 1.7 {
443     struct source *bestbuddy = NULL;
444     double d, r, mindist = MAXBUDDY;
445    
446     r = sqrt(s->dom/PI);
447     for ( ; l != NULL; l = l->next) {
448     d = sqrt(dist2(l->dir, s->dir)) - sqrt(l->dom/PI) - r;
449     if (d < mindist) {
450     bestbuddy = l;
451     mindist = d;
452     }
453     }
454     return(bestbuddy);
455     }
456    
457    
458 schorsch 2.5 extern void
459     absorb_specks(void) /* eliminate too-small sources */
460 greg 1.7 {
461     struct source head, *buddy;
462     register struct source *last, *this;
463    
464     if (verbose)
465     fprintf(stderr, "%s: absorbing small sources...\n", progname);
466     head.next = donelist;
467     last = &head;
468 greg 1.11 for (this = head.next; this != NULL; this = this->next)
469 greg 1.7 if (TOOSMALL(this)) {
470     last->next = this->next;
471 greg 1.11 buddy = findbuddy(this, head.next);
472 greg 1.7 if (buddy != NULL)
473     mergesource(buddy, this);
474     else
475     absorb(this);
476     this = last;
477     } else
478     last = this;
479     donelist = head.next;
480     }
481    
482    
483 schorsch 2.5 static void
484     absorb( /* absorb a source into indirect */
485     register struct source *s
486     )
487 greg 1.7 {
488     FVECT dir;
489 greg 1.10 double d;
490     register int i;
491 greg 1.7
492     for (i = 0; i < nglardirs; i++) {
493     spinvector(dir, ourview.vdir, ourview.vup, indirect[i].theta);
494 greg 1.10 d = DOT(dir,s->dir)*s->dom*(sampdens*sampdens);
495     if (d <= 0.0)
496 greg 1.7 continue;
497 greg 1.10 indirect[i].sum += d * s->brt;
498     indirect[i].n += d;
499 greg 1.7 }
500 greg 1.15 freespans(s);
501 greg 2.4 free((void *)s);
502 greg 1.15 }
503    
504    
505 schorsch 2.5 static void
506     freespans( /* free spans associated with source */
507     struct source *sp
508     )
509 greg 1.15 {
510     register struct srcspan *ss;
511    
512 greg 2.3 while ((ss = sp->first) != NULL) {
513     sp->first = ss->next;
514 greg 2.4 free((void *)ss);
515 greg 2.3 }
516 greg 1.1 }