ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/util/glaresrc.c
Revision: 2.7
Committed: Tue Nov 29 20:45:21 2022 UTC (15 months, 4 weeks ago) by greg
Content type: text/plain
Branch: MAIN
CVS Tags: rad5R4, HEAD
Changes since 2.6: +26 -27 lines
Log Message:
refactor(findglare): Removed outdated type modifiers

File Contents

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