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

# Content
1 #ifndef lint
2 static const char RCSid[] = "$Id: glaresrc.c,v 2.5 2004/01/02 12:51:54 schorsch Exp $";
3 #endif
4 /*
5 * Gather samples and compute glare sources.
6 */
7
8 #include "glare.h"
9 #include "linregr.h"
10
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 void pict_stats(void);
18
19 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 {
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 extern void
59 analyze(void) /* analyze our scene */
60 {
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 close_sources(v);
70 #ifndef DEBUG
71 if (verbose) {
72 fprintf(stderr, "%s: analyzing... %3ld%%\r",
73 progname, 100L*(vsize-v)/(2*vsize));
74 fflush(stderr);
75 }
76 #endif
77 getviewspan(v, spanbr);
78 left = hsize + 1;
79 for (h = -hsize; h <= hsize; h++) {
80 if (spanbr[h+hsize] < 0.0) { /* off view */
81 if (left < h) {
82 addsrcspan(newspan(left,h,v,spanbr));
83 left = hsize + 1;
84 }
85 continue;
86 }
87 if (spanbr[h+hsize] > threshold) { /* in source */
88 if (left > h)
89 left = h;
90 } else { /* out of source */
91 if (left < h) {
92 addsrcspan(newspan(left,h,v,spanbr));
93 left = hsize + 1;
94 }
95 addindirect(h, v, spanbr[h+hsize]);
96 }
97 }
98 if (left < h)
99 addsrcspan(newspan(left,h,v,spanbr));
100 }
101 free((void *)spanbr);
102 close_allsrcs();
103 }
104
105
106 static void
107 addindirect( /* add brightness to indirect illuminances */
108 int h,
109 int v,
110 double br
111 )
112 {
113 double tanb, d;
114 int hl;
115 register int i;
116
117 hl = hlim(v);
118 if (h <= -hl) { /* left region */
119 d = (double)(-h-hl)/sampdens;
120 if (d >= 1.0-FTINY)
121 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 indirect[i].n += d;
128 }
129 }
130 return;
131 }
132 if (h >= hl) { /* right region */
133 d = (double)(-h+hl)/sampdens;
134 if (d <= -1.0+FTINY)
135 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 indirect[i].n += d;
142 }
143 }
144 return;
145 }
146 /* central region */
147 for (i = 0; i < nglardirs; i++) {
148 d = cos(h_theta(h,v) - indirect[i].theta);
149 if (d > 0.0) {
150 indirect[i].sum += d * br;
151 indirect[i].n += d;
152 }
153 }
154 }
155
156
157 extern void
158 comp_thresh(void) /* compute glare threshold */
159 {
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 for (v = vsize; v >= -vsize; v -= TSAMPSTEP) {
170 for (h = -hsize; h <= hsize; h += TSAMPSTEP) {
171 if ((br = getviewpix(h, v)) < 0.0)
172 continue;
173 brsum += br;
174 nsamps++;
175 }
176 }
177 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 if (verbose) {
187 #ifdef DEBUG
188 pict_stats();
189 #endif
190 fprintf(stderr,
191 "%s: threshold set to %f cd/m2 from %d samples\n",
192 progname, threshold, nsamps);
193 }
194 }
195
196
197 static void
198 addsrcspan( /* add new source span to our list */
199 struct srcspan *nss
200 )
201 {
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 last->next = this->next;
215 mergesource(cs, this);
216 this = last;
217 }
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 cs->dom = 0.0;
228 cs->first = NULL;
229 cs->next = curlist;
230 curlist = cs;
231 }
232 nss->next = cs->first;
233 cs->first = nss;
234 }
235
236
237 static void
238 mergesource( /* merge source ap into source sp */
239 struct source *sp,
240 struct source *ap
241 )
242 {
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 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 free((void *)ap);
271 }
272
273
274 static void
275 close_sources( /* close sources above v */
276 int v
277 )
278 {
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 static void
296 close_allsrcs(void) /* done with everything */
297 {
298 register struct source *this, *next;
299
300 this = curlist;
301 while (this != NULL) {
302 next = this->next;
303 donesource(this);
304 this = next;
305 }
306 curlist = NULL;
307 }
308
309
310 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 {
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 static struct source *
342 splitsource( /* divide source in two if it's big and long */
343 struct source *so
344 )
345 {
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 static void
395 donesource( /* finished with this source */
396 register struct source *sp
397 )
398 {
399 struct source *newsrc;
400 register struct srcspan *ss;
401 int h, n;
402 double hsum, vsum, d;
403
404 while ((newsrc = splitsource(sp)) != NULL) /* split it? */
405 donesource(newsrc);
406 sp->dom = 0.0;
407 hsum = vsum = 0.0;
408 sp->brt = 0.0;
409 n = 0;
410 for (ss = sp->first; ss != NULL; ss = ss->next) {
411 sp->brt += ss->brsum;
412 n += ss->r - ss->l;
413 for (h = ss->l; h < ss->r; h++) {
414 d = pixsize(h, ss->v);
415 hsum += d*h;
416 vsum += d*ss->v;
417 sp->dom += d;
418 }
419 }
420 freespans(sp);
421 if (sp->dom <= FTINY) { /* must be right at edge of image */
422 free((void *)sp);
423 return;
424 }
425 sp->brt /= (double)n;
426 compdir(sp->dir, (int)(hsum/sp->dom), (int)(vsum/sp->dom));
427 sp->next = donelist;
428 donelist = sp;
429 if (verbose)
430 fprintf(stderr,
431 "%s: source at [%.3f,%.3f,%.3f], dw %.5f, br %.1f (%d samps)\n",
432 progname, sp->dir[0], sp->dir[1], sp->dir[2],
433 sp->dom, sp->brt, n);
434 }
435
436
437 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 {
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 extern void
459 absorb_specks(void) /* eliminate too-small sources */
460 {
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 for (this = head.next; this != NULL; this = this->next)
469 if (TOOSMALL(this)) {
470 last->next = this->next;
471 buddy = findbuddy(this, head.next);
472 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 static void
484 absorb( /* absorb a source into indirect */
485 register struct source *s
486 )
487 {
488 FVECT dir;
489 double d;
490 register int i;
491
492 for (i = 0; i < nglardirs; i++) {
493 spinvector(dir, ourview.vdir, ourview.vup, indirect[i].theta);
494 d = DOT(dir,s->dir)*s->dom*(sampdens*sampdens);
495 if (d <= 0.0)
496 continue;
497 indirect[i].sum += d * s->brt;
498 indirect[i].n += d;
499 }
500 freespans(s);
501 free((void *)s);
502 }
503
504
505 static void
506 freespans( /* free spans associated with source */
507 struct source *sp
508 )
509 {
510 register struct srcspan *ss;
511
512 while ((ss = sp->first) != NULL) {
513 sp->first = ss->next;
514 free((void *)ss);
515 }
516 }