35 |
|
RAY *rp; /* originating ray sample */ |
36 |
|
int ns; /* number of samples per axis */ |
37 |
|
int sampOK; /* acquired full sample set? */ |
38 |
+ |
int atyp; /* RAMBIENT or TAMBIENT */ |
39 |
|
SCOLOR acoef; /* division contribution coefficient */ |
40 |
|
SCOLOR acol; /* accumulated color */ |
41 |
+ |
FVECT onrm; /* oriented unperturbed surface normal */ |
42 |
|
FVECT ux, uy; /* tangent axis unit vectors */ |
43 |
|
AMBSAMP sa[1]; /* sample array (extends struct) */ |
44 |
|
} AMBHEMI; /* ambient sample hemisphere */ |
108 |
|
setscolor(ar.rcoef, AVGREFL, AVGREFL, AVGREFL); |
109 |
|
else |
110 |
|
copyscolor(ar.rcoef, hp->acoef); |
111 |
< |
if (rayorigin(&ar, AMBIENT, hp->rp, ar.rcoef) < 0) |
111 |
> |
if (rayorigin(&ar, hp->atyp, hp->rp, ar.rcoef) < 0) |
112 |
|
return(0); |
113 |
|
if (ambacc > FTINY) { |
114 |
|
smultscolor(ar.rcoef, hp->acoef); |
124 |
|
for (ii = 3; ii--; ) |
125 |
|
ar.rdir[ii] = spt[0]*hp->ux[ii] + |
126 |
|
spt[1]*hp->uy[ii] + |
127 |
< |
zd*hp->rp->ron[ii]; |
127 |
> |
zd*hp->onrm[ii]; |
128 |
|
checknorm(ar.rdir); |
129 |
|
/* avoid coincident samples */ |
130 |
|
if (!n && ambcollision(hp, i, j, ar.rdir)) { |
162 |
|
static float * |
163 |
|
getambdiffs(AMBHEMI *hp) |
164 |
|
{ |
165 |
< |
const double normf = 1./bright(hp->acoef); |
165 |
> |
const double normf = 1./(pbright(hp->acoef) + FTINY); |
166 |
|
float *earr = (float *)calloc(hp->ns*hp->ns, sizeof(float)); |
167 |
< |
float *ep; |
167 |
> |
float *ep, *earr2; |
168 |
|
AMBSAMP *ap; |
169 |
|
double b, b1, d2; |
170 |
|
int i, j; |
198 |
|
ep[-hp->ns-1] += d2; |
199 |
|
} |
200 |
|
/* correct for number of neighbors */ |
201 |
< |
earr[0] *= 8./3.; |
202 |
< |
earr[hp->ns-1] *= 8./3.; |
203 |
< |
earr[(hp->ns-1)*hp->ns] *= 8./3.; |
204 |
< |
earr[(hp->ns-1)*hp->ns + hp->ns-1] *= 8./3.; |
201 |
> |
earr[0] *= 6./3.; |
202 |
> |
earr[hp->ns-1] *= 6./3.; |
203 |
> |
earr[(hp->ns-1)*hp->ns] *= 6./3.; |
204 |
> |
earr[(hp->ns-1)*hp->ns + hp->ns-1] *= 6./3.; |
205 |
|
for (i = 1; i < hp->ns-1; i++) { |
206 |
< |
earr[i*hp->ns] *= 8./5.; |
207 |
< |
earr[i*hp->ns + hp->ns-1] *= 8./5.; |
206 |
> |
earr[i*hp->ns] *= 6./5.; |
207 |
> |
earr[i*hp->ns + hp->ns-1] *= 6./5.; |
208 |
|
} |
209 |
|
for (j = 1; j < hp->ns-1; j++) { |
210 |
< |
earr[j] *= 8./5.; |
211 |
< |
earr[(hp->ns-1)*hp->ns + j] *= 8./5.; |
210 |
> |
earr[j] *= 6./5.; |
211 |
> |
earr[(hp->ns-1)*hp->ns + j] *= 6./5.; |
212 |
|
} |
213 |
+ |
/* preen map to avoid cliffs */ |
214 |
+ |
earr2 = (float *)malloc(hp->ns*hp->ns*sizeof(float)); |
215 |
+ |
if (earr2 == NULL) |
216 |
+ |
return(earr); |
217 |
+ |
memcpy(earr2, earr, hp->ns*hp->ns*sizeof(float)); |
218 |
+ |
for (i = 0; i < hp->ns-1; i++) { |
219 |
+ |
float *ep2 = earr2 + i*hp->ns; |
220 |
+ |
ep = earr + i*hp->ns; |
221 |
+ |
for (j = 0; j < hp->ns-1; j++, ep2++, ep++) { |
222 |
+ |
if (ep2[1] < .5*ep2[0]) { |
223 |
+ |
ep[0] -= .125*ep2[0]; |
224 |
+ |
ep[1] += .125*ep2[0]; |
225 |
+ |
} else if (ep2[1] > 2.*ep2[0]) { |
226 |
+ |
ep[1] -= .125*ep2[1]; |
227 |
+ |
ep[0] += .125*ep2[1]; |
228 |
+ |
} |
229 |
+ |
if (ep2[hp->ns] < .5*ep2[0]) { |
230 |
+ |
ep[0] -= .125*ep2[0]; |
231 |
+ |
ep[hp->ns] += .125*ep2[0]; |
232 |
+ |
} else if (ep2[hp->ns] > 2.*ep2[0]) { |
233 |
+ |
ep[hp->ns] -= .125*ep2[hp->ns]; |
234 |
+ |
ep[0] += .125*ep2[hp->ns]; |
235 |
+ |
} |
236 |
+ |
} |
237 |
+ |
} |
238 |
+ |
free(earr2); |
239 |
|
return(earr); |
240 |
|
} |
241 |
|
|
276 |
|
double wt |
277 |
|
) |
278 |
|
{ |
279 |
+ |
int backside = (wt < 0); |
280 |
|
AMBHEMI *hp; |
281 |
|
double d; |
282 |
|
int n, i, j; |
285 |
|
if (d <= FTINY) |
286 |
|
return(NULL); |
287 |
|
/* set number of divisions */ |
288 |
+ |
if (backside) wt = -wt; |
289 |
|
if (ambacc <= FTINY && |
290 |
|
wt > (d *= 0.8*r->rweight/(ambdiv*minweight))) |
291 |
|
wt = d; /* avoid ray termination */ |
297 |
|
hp = (AMBHEMI *)malloc(sizeof(AMBHEMI) + sizeof(AMBSAMP)*(n*n - 1)); |
298 |
|
if (hp == NULL) |
299 |
|
error(SYSTEM, "out of memory in samp_hemi"); |
300 |
+ |
|
301 |
+ |
if (backside) { |
302 |
+ |
hp->atyp = TAMBIENT; |
303 |
+ |
hp->onrm[0] = -r->ron[0]; |
304 |
+ |
hp->onrm[1] = -r->ron[1]; |
305 |
+ |
hp->onrm[2] = -r->ron[2]; |
306 |
+ |
} else { |
307 |
+ |
hp->atyp = RAMBIENT; |
308 |
+ |
VCOPY(hp->onrm, r->ron); |
309 |
+ |
} |
310 |
|
hp->rp = r; |
311 |
|
hp->ns = n; |
312 |
|
scolorblack(hp->acol); |
317 |
|
d = 1.0/(n*n); |
318 |
|
scalescolor(hp->acoef, d); |
319 |
|
/* make tangent plane axes */ |
320 |
< |
if (!getperpendicular(hp->ux, r->ron, 1)) |
320 |
> |
if (!getperpendicular(hp->ux, hp->onrm, 1)) |
321 |
|
error(CONSISTENCY, "bad ray direction in samp_hemi"); |
322 |
< |
VCROSS(hp->uy, r->ron, hp->ux); |
322 |
> |
VCROSS(hp->uy, hp->onrm, hp->ux); |
323 |
|
/* sample divisions */ |
324 |
|
for (i = hp->ns; i--; ) |
325 |
|
for (j = hp->ns; j--; ) |
584 |
|
for (j = 0; j < hp->ns-1; j++) { |
585 |
|
comp_fftri(&fftr, hp, AI(hp,0,j), AI(hp,0,j+1)); |
586 |
|
if (hessrow != NULL) |
587 |
< |
comp_hessian(hessrow[j], &fftr, hp->rp->ron); |
587 |
> |
comp_hessian(hessrow[j], &fftr, hp->onrm); |
588 |
|
if (gradrow != NULL) |
589 |
< |
comp_gradient(gradrow[j], &fftr, hp->rp->ron); |
589 |
> |
comp_gradient(gradrow[j], &fftr, hp->onrm); |
590 |
|
} |
591 |
|
/* sum each row of triangles */ |
592 |
|
for (i = 0; i < hp->ns-1; i++) { |
594 |
|
FVECT gradcol; |
595 |
|
comp_fftri(&fftr, hp, AI(hp,i,0), AI(hp,i+1,0)); |
596 |
|
if (hessrow != NULL) |
597 |
< |
comp_hessian(hesscol, &fftr, hp->rp->ron); |
597 |
> |
comp_hessian(hesscol, &fftr, hp->onrm); |
598 |
|
if (gradrow != NULL) |
599 |
< |
comp_gradient(gradcol, &fftr, hp->rp->ron); |
599 |
> |
comp_gradient(gradcol, &fftr, hp->onrm); |
600 |
|
for (j = 0; j < hp->ns-1; j++) { |
601 |
|
FVECT hessdia[3]; /* compute triangle contributions */ |
602 |
|
FVECT graddia; |
606 |
|
/* diagonal (inner) edge */ |
607 |
|
comp_fftri(&fftr, hp, AI(hp,i,j+1), AI(hp,i+1,j)); |
608 |
|
if (hessrow != NULL) { |
609 |
< |
comp_hessian(hessdia, &fftr, hp->rp->ron); |
609 |
> |
comp_hessian(hessdia, &fftr, hp->onrm); |
610 |
|
rev_hessian(hesscol); |
611 |
|
add2hessian(hessian, hessrow[j], hessdia, hesscol, backg); |
612 |
|
} |
613 |
|
if (gradrow != NULL) { |
614 |
< |
comp_gradient(graddia, &fftr, hp->rp->ron); |
614 |
> |
comp_gradient(graddia, &fftr, hp->onrm); |
615 |
|
rev_gradient(gradcol); |
616 |
|
add2gradient(gradient, gradrow[j], graddia, gradcol, backg); |
617 |
|
} |
618 |
|
/* initialize edge in next row */ |
619 |
|
comp_fftri(&fftr, hp, AI(hp,i+1,j+1), AI(hp,i+1,j)); |
620 |
|
if (hessrow != NULL) |
621 |
< |
comp_hessian(hessrow[j], &fftr, hp->rp->ron); |
621 |
> |
comp_hessian(hessrow[j], &fftr, hp->onrm); |
622 |
|
if (gradrow != NULL) |
623 |
< |
comp_gradient(gradrow[j], &fftr, hp->rp->ron); |
623 |
> |
comp_gradient(gradrow[j], &fftr, hp->onrm); |
624 |
|
/* new column edge & paired triangle */ |
625 |
|
backg = back_ambval(hp, AI(hp,i+1,j+1), |
626 |
|
AI(hp,i+1,j), AI(hp,i,j+1)); |
627 |
|
comp_fftri(&fftr, hp, AI(hp,i,j+1), AI(hp,i+1,j+1)); |
628 |
|
if (hessrow != NULL) { |
629 |
< |
comp_hessian(hesscol, &fftr, hp->rp->ron); |
629 |
> |
comp_hessian(hesscol, &fftr, hp->onrm); |
630 |
|
rev_hessian(hessdia); |
631 |
|
add2hessian(hessian, hessrow[j], hessdia, hesscol, backg); |
632 |
|
if (i < hp->ns-2) |
633 |
|
rev_hessian(hessrow[j]); |
634 |
|
} |
635 |
|
if (gradrow != NULL) { |
636 |
< |
comp_gradient(gradcol, &fftr, hp->rp->ron); |
636 |
> |
comp_gradient(gradcol, &fftr, hp->onrm); |
637 |
|
rev_gradient(graddia); |
638 |
|
add2gradient(gradient, gradrow[j], graddia, gradcol, backg); |
639 |
|
if (i < hp->ns-2) |
669 |
|
/* use vector for azimuth + 90deg */ |
670 |
|
VSUB(vd, ap->p, hp->rp->rop); |
671 |
|
/* brightness over cosine factor */ |
672 |
< |
gfact = ap->v[0] / DOT(hp->rp->ron, vd); |
672 |
> |
gfact = ap->v[0] / DOT(hp->onrm, vd); |
673 |
|
/* sine = proj_radius/vd_length */ |
674 |
|
dgsum[0] -= DOT(uv[1], vd) * gfact; |
675 |
|
dgsum[1] += DOT(uv[0], vd) * gfact; |
727 |
|
doambient( /* compute ambient component */ |
728 |
|
SCOLOR rcol, /* input/output color */ |
729 |
|
RAY *r, |
730 |
< |
double wt, |
730 |
> |
double wt, /* negative for back side */ |
731 |
|
FVECT uv[2], /* returned (optional) */ |
732 |
|
float ra[2], /* returned (optional) */ |
733 |
|
float pg[2], /* returned (optional) */ |
759 |
|
free(hp); /* Hessian not requested/possible */ |
760 |
|
return(-1); /* value-only return value */ |
761 |
|
} |
762 |
< |
if ((d = scolor_photopic(rcol)) > FTINY) { |
763 |
< |
d = 0.99*(hp->ns*hp->ns)/d; /* normalize Y values */ |
762 |
> |
if ((d = scolor_mean(rcol)) > FTINY) { |
763 |
> |
d = 0.99*(hp->ns*hp->ns)/d; /* normalize avg. values */ |
764 |
|
K = 0.01; |
765 |
|
} else { /* or fall back on geometric Hessian */ |
766 |
|
K = 1.0; |
794 |
|
if (ra[1] < minarad) |
795 |
|
ra[1] = minarad; |
796 |
|
} |
797 |
< |
ra[0] *= d = 1.0/sqrt(wt); |
797 |
> |
ra[0] *= d = 1.0/sqrt(fabs(wt)); |
798 |
|
if ((ra[1] *= d) > 2.0*ra[0]) |
799 |
|
ra[1] = 2.0*ra[0]; |
800 |
|
if (ra[1] > maxarad) { |