| 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 |
#include "linregr.h"
|
| 13 |
|
| 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 |
close_sources(v);
|
| 53 |
#ifndef DEBUG
|
| 54 |
if (verbose)
|
| 55 |
fprintf(stderr, "%s: analyzing... %3ld%%\r",
|
| 56 |
progname, 100L*(vsize-v)/(2*vsize));
|
| 57 |
#endif
|
| 58 |
getviewspan(v, spanbr);
|
| 59 |
left = hsize + 1;
|
| 60 |
for (h = -hsize; h <= hsize; h++) {
|
| 61 |
if (spanbr[h+hsize] < 0.0) { /* off view */
|
| 62 |
if (left < h) {
|
| 63 |
addsrcspan(newspan(left,h,v,spanbr));
|
| 64 |
left = hsize + 1;
|
| 65 |
}
|
| 66 |
continue;
|
| 67 |
}
|
| 68 |
if (spanbr[h+hsize] > threshold) { /* in source */
|
| 69 |
if (left > h)
|
| 70 |
left = h;
|
| 71 |
} else { /* out of source */
|
| 72 |
if (left < h) {
|
| 73 |
addsrcspan(newspan(left,h,v,spanbr));
|
| 74 |
left = hsize + 1;
|
| 75 |
}
|
| 76 |
addindirect(h, v, spanbr[h+hsize]);
|
| 77 |
}
|
| 78 |
}
|
| 79 |
if (left < h)
|
| 80 |
addsrcspan(newspan(left,h,v,spanbr));
|
| 81 |
}
|
| 82 |
free((char *)spanbr);
|
| 83 |
close_allsrcs();
|
| 84 |
}
|
| 85 |
|
| 86 |
|
| 87 |
addindirect(h, v, br) /* add brightness to indirect illuminances */
|
| 88 |
int h, v;
|
| 89 |
double br;
|
| 90 |
{
|
| 91 |
double tanb, d;
|
| 92 |
int hl;
|
| 93 |
register int i;
|
| 94 |
|
| 95 |
hl = hlim(v);
|
| 96 |
if (h <= -hl) { /* left region */
|
| 97 |
d = (double)(-h-hl)/sampdens;
|
| 98 |
if (d >= 1.0-FTINY)
|
| 99 |
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 |
indirect[i].n += d;
|
| 106 |
}
|
| 107 |
}
|
| 108 |
return;
|
| 109 |
}
|
| 110 |
if (h >= hl) { /* right region */
|
| 111 |
d = (double)(-h+hl)/sampdens;
|
| 112 |
if (d <= -1.0+FTINY)
|
| 113 |
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 |
indirect[i].n += d;
|
| 120 |
}
|
| 121 |
}
|
| 122 |
return;
|
| 123 |
}
|
| 124 |
/* central region */
|
| 125 |
for (i = 0; i < nglardirs; i++) {
|
| 126 |
d = cos(h_theta(h,v) - indirect[i].theta);
|
| 127 |
if (d > 0.0) {
|
| 128 |
indirect[i].sum += d * br;
|
| 129 |
indirect[i].n += d;
|
| 130 |
}
|
| 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 |
for (v = vsize; v >= -vsize; v -= TSAMPSTEP) {
|
| 147 |
for (h = -hsize; h <= hsize; h += TSAMPSTEP) {
|
| 148 |
if ((br = getviewpix(h, v)) < 0.0)
|
| 149 |
continue;
|
| 150 |
brsum += br;
|
| 151 |
nsamps++;
|
| 152 |
}
|
| 153 |
}
|
| 154 |
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 |
if (verbose) {
|
| 164 |
#ifdef DEBUG
|
| 165 |
pict_stats();
|
| 166 |
#endif
|
| 167 |
fprintf(stderr,
|
| 168 |
"%s: threshold set to %f cd/m2 from %d samples\n",
|
| 169 |
progname, threshold, nsamps);
|
| 170 |
}
|
| 171 |
}
|
| 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 |
last->next = this->next;
|
| 190 |
mergesource(cs, this);
|
| 191 |
this = last;
|
| 192 |
}
|
| 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 |
cs->dom = 0.0;
|
| 203 |
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 |
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 |
}
|
| 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 |
register struct source *this, *next;
|
| 268 |
|
| 269 |
this = curlist;
|
| 270 |
while (this != NULL) {
|
| 271 |
next = this->next;
|
| 272 |
donesource(this);
|
| 273 |
this = next;
|
| 274 |
}
|
| 275 |
curlist = NULL;
|
| 276 |
}
|
| 277 |
|
| 278 |
|
| 279 |
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 |
donesource(sp) /* finished with this source */
|
| 360 |
register struct source *sp;
|
| 361 |
{
|
| 362 |
struct source *newsrc;
|
| 363 |
register struct srcspan *ss;
|
| 364 |
int h, n;
|
| 365 |
double hsum, vsum, d;
|
| 366 |
|
| 367 |
while ((newsrc = splitsource(sp)) != NULL) /* split it? */
|
| 368 |
donesource(newsrc);
|
| 369 |
sp->dom = 0.0;
|
| 370 |
hsum = vsum = 0.0;
|
| 371 |
sp->brt = 0.0;
|
| 372 |
n = 0;
|
| 373 |
for (ss = sp->first; ss != NULL; ss = ss->next) {
|
| 374 |
sp->brt += ss->brsum;
|
| 375 |
n += ss->r - ss->l;
|
| 376 |
for (h = ss->l; h < ss->r; h++) {
|
| 377 |
d = pixsize(h, ss->v);
|
| 378 |
hsum += d*h;
|
| 379 |
vsum += d*ss->v;
|
| 380 |
sp->dom += d;
|
| 381 |
}
|
| 382 |
}
|
| 383 |
sp->brt /= (double)n;
|
| 384 |
compdir(sp->dir, (int)(hsum/sp->dom), (int)(vsum/sp->dom));
|
| 385 |
freespans(sp);
|
| 386 |
sp->next = donelist;
|
| 387 |
donelist = sp;
|
| 388 |
if (verbose)
|
| 389 |
fprintf(stderr,
|
| 390 |
"%s: source at [%.3f,%.3f,%.3f], dw %.5f, br %.1f (%d samps)\n",
|
| 391 |
progname, sp->dir[0], sp->dir[1], sp->dir[2],
|
| 392 |
sp->dom, sp->brt, n);
|
| 393 |
}
|
| 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 |
for (this = head.next; this != NULL; this = this->next)
|
| 425 |
if (TOOSMALL(this)) {
|
| 426 |
last->next = this->next;
|
| 427 |
buddy = findbuddy(this, head.next);
|
| 428 |
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 |
double d;
|
| 444 |
register int i;
|
| 445 |
|
| 446 |
for (i = 0; i < nglardirs; i++) {
|
| 447 |
spinvector(dir, ourview.vdir, ourview.vup, indirect[i].theta);
|
| 448 |
d = DOT(dir,s->dir)*s->dom*(sampdens*sampdens);
|
| 449 |
if (d <= 0.0)
|
| 450 |
continue;
|
| 451 |
indirect[i].sum += d * s->brt;
|
| 452 |
indirect[i].n += d;
|
| 453 |
}
|
| 454 |
freespans(s);
|
| 455 |
free((char *)s);
|
| 456 |
}
|
| 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 |
}
|
