ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/util/glaresrc.c
Revision: 1.15
Committed: Tue Apr 30 16:39:39 1991 UTC (32 years, 11 months ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 1.14: +97 -4 lines
Log Message:
added routines for splitting long sources

File Contents

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