24 |
|
#ifndef MINADIV |
25 |
|
#define MINADIV 7 /* minimum # divisions in each dimension */ |
26 |
|
#endif |
27 |
+ |
#ifndef MINSDIST |
28 |
+ |
#define MINSDIST 0.25 /* def. min. spacing = 1/4th division */ |
29 |
+ |
#endif |
30 |
|
|
28 |
– |
extern void SDsquare2disk(double ds[2], double seedx, double seedy); |
29 |
– |
|
31 |
|
typedef struct { |
31 |
– |
COLOR v; /* hemisphere sample value */ |
32 |
– |
float d; /* reciprocal distance */ |
32 |
|
FVECT p; /* intersection point */ |
33 |
+ |
float d; /* reciprocal distance */ |
34 |
+ |
SCOLOR v; /* hemisphere sample value */ |
35 |
|
} AMBSAMP; /* sample value */ |
36 |
|
|
37 |
|
typedef struct { |
38 |
|
RAY *rp; /* originating ray sample */ |
39 |
|
int ns; /* number of samples per axis */ |
40 |
|
int sampOK; /* acquired full sample set? */ |
41 |
< |
COLOR acoef; /* division contribution coefficient */ |
42 |
< |
double acol[3]; /* accumulated color */ |
41 |
> |
int atyp; /* RAMBIENT or TAMBIENT */ |
42 |
> |
SCOLOR acoef; /* division contribution coefficient */ |
43 |
> |
SCOLOR acol; /* accumulated color */ |
44 |
> |
FVECT onrm; /* oriented unperturbed surface normal */ |
45 |
|
FVECT ux, uy; /* tangent axis unit vectors */ |
46 |
|
AMBSAMP sa[1]; /* sample array (extends struct) */ |
47 |
|
} AMBHEMI; /* ambient sample hemisphere */ |
65 |
|
{ |
66 |
|
double cos_thresh; |
67 |
|
int ii, jj; |
68 |
< |
/* min. spacing = 1/4th division */ |
69 |
< |
cos_thresh = (PI/4.)/(double)hp->ns; |
68 |
> |
|
69 |
> |
cos_thresh = (PI*MINSDIST)/(double)hp->ns; |
70 |
|
cos_thresh = 1. - .5*cos_thresh*cos_thresh; |
71 |
|
/* check existing neighbors */ |
72 |
|
for (ii = i-1; ii <= i+1; ii++) { |
92 |
|
} |
93 |
|
|
94 |
|
|
95 |
+ |
#define XLOTSIZ 251 /* size of used car lot */ |
96 |
+ |
#define CFIRST 0 /* first corner */ |
97 |
+ |
#define COTHER (CFIRST+4) /* non-corner sample */ |
98 |
+ |
#define CMAXTARGET (int)(XLOTSIZ*MINSDIST/(1-MINSDIST)) |
99 |
+ |
#define CXCOPY(d,s) (excharr[d][0]=excharr[s][0], excharr[d][1]=excharr[s][1]) |
100 |
+ |
|
101 |
|
static int |
102 |
+ |
psample_class(double ss[2]) /* classify patch sample */ |
103 |
+ |
{ |
104 |
+ |
if (ss[0] < MINSDIST) { |
105 |
+ |
if (ss[1] < MINSDIST) |
106 |
+ |
return(CFIRST); |
107 |
+ |
if (ss[1] > 1.-MINSDIST) |
108 |
+ |
return(CFIRST+2); |
109 |
+ |
} else if (ss[0] > 1.-MINSDIST) { |
110 |
+ |
if (ss[1] < MINSDIST) |
111 |
+ |
return(CFIRST+1); |
112 |
+ |
if (ss[1] > 1.-MINSDIST) |
113 |
+ |
return(CFIRST+3); |
114 |
+ |
} |
115 |
+ |
return(COTHER); /* not in a corner */ |
116 |
+ |
} |
117 |
+ |
|
118 |
+ |
static void |
119 |
+ |
trade_patchsamp(double ss[2]) /* trade in problem patch position */ |
120 |
+ |
{ |
121 |
+ |
static float excharr[XLOTSIZ][2]; |
122 |
+ |
static short gterm[COTHER+1]; |
123 |
+ |
double srep[2]; |
124 |
+ |
int sclass, rclass; |
125 |
+ |
int x; |
126 |
+ |
/* reset on corner overload */ |
127 |
+ |
if (gterm[COTHER-1] >= (CMAXTARGET+XLOTSIZ)/2) |
128 |
+ |
memset(gterm, 0, sizeof(gterm)); |
129 |
+ |
/* (re-)initialize? */ |
130 |
+ |
while (gterm[COTHER] < XLOTSIZ) { |
131 |
+ |
excharr[gterm[COTHER]][0] = frandom(); |
132 |
+ |
excharr[gterm[COTHER]][1] = frandom(); |
133 |
+ |
++gterm[COTHER]; |
134 |
+ |
} /* get trade-in candidate... */ |
135 |
+ |
sclass = psample_class(ss); /* submitted corner or not? */ |
136 |
+ |
switch (sclass) { |
137 |
+ |
case COTHER: /* trade mid-edge with corner/any */ |
138 |
+ |
x = irandom( gterm[COTHER-1] > CMAXTARGET |
139 |
+ |
? gterm[COTHER-1] : XLOTSIZ ); |
140 |
+ |
break; |
141 |
+ |
case CFIRST: /* kick out of first corner */ |
142 |
+ |
x = gterm[CFIRST] + irandom(XLOTSIZ - gterm[CFIRST]); |
143 |
+ |
break; |
144 |
+ |
default: /* kick out of 2nd-4th corner */ |
145 |
+ |
x = irandom(XLOTSIZ - (gterm[sclass] - gterm[sclass-1])); |
146 |
+ |
x += (x >= gterm[sclass-1])*(gterm[sclass] - gterm[sclass-1]); |
147 |
+ |
break; |
148 |
+ |
} |
149 |
+ |
srep[0] = excharr[x][0]; /* save selected trade output */ |
150 |
+ |
srep[1] = excharr[x][1]; |
151 |
+ |
/* adjust our lot groups */ |
152 |
+ |
for (rclass = CFIRST; rclass < COTHER; rclass++) |
153 |
+ |
if (x < gterm[rclass]) |
154 |
+ |
break; |
155 |
+ |
if (sclass < rclass) { /* submitted group before replacement? */ |
156 |
+ |
CXCOPY(x, gterm[rclass-1]); |
157 |
+ |
while (--rclass > sclass) { |
158 |
+ |
CXCOPY(gterm[rclass], gterm[rclass-1]); |
159 |
+ |
++gterm[rclass]; |
160 |
+ |
} |
161 |
+ |
x = gterm[sclass]++; |
162 |
+ |
} else if (sclass > rclass) { /* submitted group after replacement? */ |
163 |
+ |
--gterm[rclass]; |
164 |
+ |
CXCOPY(x, gterm[rclass]); |
165 |
+ |
while (++rclass < sclass) { |
166 |
+ |
--gterm[rclass]; |
167 |
+ |
CXCOPY(gterm[rclass-1], gterm[rclass]); |
168 |
+ |
} |
169 |
+ |
x = gterm[sclass-1]; |
170 |
+ |
} |
171 |
+ |
excharr[x][0] = ss[0]; /* complete the transaction */ |
172 |
+ |
excharr[x][1] = ss[1]; |
173 |
+ |
ss[0] = srep[0]; |
174 |
+ |
ss[1] = srep[1]; |
175 |
+ |
} |
176 |
+ |
|
177 |
+ |
#undef CXCOPY |
178 |
+ |
#undef XLOTSIZ |
179 |
+ |
#undef COTHER |
180 |
+ |
#undef CFIRST |
181 |
+ |
|
182 |
+ |
|
183 |
+ |
static int |
184 |
|
ambsample( /* initial ambient division sample */ |
185 |
|
AMBHEMI *hp, |
186 |
|
int i, |
191 |
|
AMBSAMP *ap = &ambsam(hp,i,j); |
192 |
|
RAY ar; |
193 |
|
int hlist[3], ii; |
194 |
< |
double spt[2], zd; |
194 |
> |
double ss[2]; |
195 |
> |
RREAL spt[2]; |
196 |
> |
double zd; |
197 |
|
/* generate hemispherical sample */ |
198 |
|
/* ambient coefficient for weight */ |
199 |
|
if (ambacc > FTINY) |
200 |
< |
setcolor(ar.rcoef, AVGREFL, AVGREFL, AVGREFL); |
200 |
> |
setscolor(ar.rcoef, AVGREFL, AVGREFL, AVGREFL); |
201 |
|
else |
202 |
< |
copycolor(ar.rcoef, hp->acoef); |
203 |
< |
if (rayorigin(&ar, AMBIENT, hp->rp, ar.rcoef) < 0) |
202 |
> |
copyscolor(ar.rcoef, hp->acoef); |
203 |
> |
if (rayorigin(&ar, hp->atyp, hp->rp, ar.rcoef) < 0) |
204 |
|
return(0); |
205 |
|
if (ambacc > FTINY) { |
206 |
< |
multcolor(ar.rcoef, hp->acoef); |
207 |
< |
scalecolor(ar.rcoef, 1./AVGREFL); |
206 |
> |
smultscolor(ar.rcoef, hp->acoef); |
207 |
> |
scalescolor(ar.rcoef, 1./AVGREFL); |
208 |
|
} |
209 |
|
hlist[0] = hp->rp->rno; |
210 |
< |
hlist[1] = j; |
211 |
< |
hlist[2] = i; |
212 |
< |
multisamp(spt, 2, urand(ilhash(hlist,3)+n)); |
213 |
< |
resample: |
214 |
< |
SDsquare2disk(spt, (j+spt[1])/hp->ns, (i+spt[0])/hp->ns); |
210 |
> |
hlist[1] = AI(hp,i,j); |
211 |
> |
hlist[2] = samplendx; |
212 |
> |
multisamp(ss, 2, urand(ilhash(hlist,3)+n)); |
213 |
> |
patch_redo: |
214 |
> |
square2disk(spt, (j+ss[1])/hp->ns, (i+ss[0])/hp->ns); |
215 |
|
zd = sqrt(1. - spt[0]*spt[0] - spt[1]*spt[1]); |
216 |
|
for (ii = 3; ii--; ) |
217 |
|
ar.rdir[ii] = spt[0]*hp->ux[ii] + |
218 |
|
spt[1]*hp->uy[ii] + |
219 |
< |
zd*hp->rp->ron[ii]; |
219 |
> |
zd*hp->onrm[ii]; |
220 |
|
checknorm(ar.rdir); |
221 |
|
/* avoid coincident samples */ |
222 |
< |
if (!n && ambcollision(hp, i, j, ar.rdir)) { |
223 |
< |
spt[0] = frandom(); spt[1] = frandom(); |
224 |
< |
goto resample; /* reject this sample */ |
222 |
> |
if (!n & (hp->ns >= 4) && ambcollision(hp, i, j, ar.rdir)) { |
223 |
> |
trade_patchsamp(ss); |
224 |
> |
goto patch_redo; |
225 |
|
} |
226 |
|
dimlist[ndims++] = AI(hp,i,j) + 90171; |
227 |
|
rayvalue(&ar); /* evaluate ray */ |
229 |
|
zd = raydistance(&ar); |
230 |
|
if (zd <= FTINY) |
231 |
|
return(0); /* should never happen */ |
232 |
< |
multcolor(ar.rcol, ar.rcoef); /* apply coefficient */ |
232 |
> |
smultscolor(ar.rcol, ar.rcoef); /* apply coefficient */ |
233 |
|
if (zd*ap->d < 1.0) /* new/closer distance? */ |
234 |
|
ap->d = 1.0/zd; |
235 |
|
if (!n) { /* record first vertex & value */ |
236 |
|
if (zd > 10.0*thescene.cusize + 1000.) |
237 |
|
zd = 10.0*thescene.cusize + 1000.; |
238 |
|
VSUM(ap->p, ar.rorg, ar.rdir, zd); |
239 |
< |
copycolor(ap->v, ar.rcol); |
239 |
> |
copyscolor(ap->v, ar.rcol); |
240 |
|
} else { /* else update recorded value */ |
241 |
< |
hp->acol[RED] -= colval(ap->v,RED); |
149 |
< |
hp->acol[GRN] -= colval(ap->v,GRN); |
150 |
< |
hp->acol[BLU] -= colval(ap->v,BLU); |
241 |
> |
sopscolor(hp->acol, -=, ap->v); |
242 |
|
zd = 1.0/(double)(n+1); |
243 |
< |
scalecolor(ar.rcol, zd); |
243 |
> |
scalescolor(ar.rcol, zd); |
244 |
|
zd *= (double)n; |
245 |
< |
scalecolor(ap->v, zd); |
246 |
< |
addcolor(ap->v, ar.rcol); |
245 |
> |
scalescolor(ap->v, zd); |
246 |
> |
saddscolor(ap->v, ar.rcol); |
247 |
|
} |
248 |
< |
addcolor(hp->acol, ap->v); /* add to our sum */ |
248 |
> |
saddscolor(hp->acol, ap->v); /* add to our sum */ |
249 |
|
return(1); |
250 |
|
} |
251 |
|
|
254 |
|
static float * |
255 |
|
getambdiffs(AMBHEMI *hp) |
256 |
|
{ |
257 |
< |
const double normf = 1./bright(hp->acoef); |
258 |
< |
float *earr = (float *)calloc(hp->ns*hp->ns, sizeof(float)); |
257 |
> |
const double normf = 1./(pbright(hp->acoef) + FTINY); |
258 |
> |
float *earr = (float *)calloc(2*hp->ns*hp->ns, sizeof(float)); |
259 |
|
float *ep; |
260 |
|
AMBSAMP *ap; |
261 |
|
double b, b1, d2; |
264 |
|
if (earr == NULL) /* out of memory? */ |
265 |
|
return(NULL); |
266 |
|
/* sum squared neighbor diffs */ |
267 |
< |
for (ap = hp->sa, ep = earr, i = 0; i < hp->ns; i++) |
267 |
> |
ap = hp->sa; |
268 |
> |
ep = earr + hp->ns*hp->ns; /* original estimates to scratch */ |
269 |
> |
for (i = 0; i < hp->ns; i++) |
270 |
|
for (j = 0; j < hp->ns; j++, ap++, ep++) { |
271 |
< |
b = bright(ap[0].v); |
271 |
> |
b = pbright(ap[0].v); |
272 |
|
if (i) { /* from above */ |
273 |
< |
b1 = bright(ap[-hp->ns].v); |
273 |
> |
b1 = pbright(ap[-hp->ns].v); |
274 |
|
d2 = b - b1; |
275 |
< |
d2 *= d2*normf/(b + b1); |
275 |
> |
d2 *= d2*normf/(b + b1 + FTINY); |
276 |
|
ep[0] += d2; |
277 |
|
ep[-hp->ns] += d2; |
278 |
|
} |
279 |
|
if (!j) continue; |
280 |
|
/* from behind */ |
281 |
< |
b1 = bright(ap[-1].v); |
281 |
> |
b1 = pbright(ap[-1].v); |
282 |
|
d2 = b - b1; |
283 |
< |
d2 *= d2*normf/(b + b1); |
283 |
> |
d2 *= d2*normf/(b + b1 + FTINY); |
284 |
|
ep[0] += d2; |
285 |
|
ep[-1] += d2; |
286 |
|
if (!i) continue; |
287 |
|
/* diagonal */ |
288 |
< |
b1 = bright(ap[-hp->ns-1].v); |
288 |
> |
b1 = pbright(ap[-hp->ns-1].v); |
289 |
|
d2 = b - b1; |
290 |
< |
d2 *= d2*normf/(b + b1); |
290 |
> |
d2 *= d2*normf/(b + b1 + FTINY); |
291 |
|
ep[0] += d2; |
292 |
|
ep[-hp->ns-1] += d2; |
293 |
|
} |
294 |
|
/* correct for number of neighbors */ |
295 |
< |
earr[0] *= 8./3.; |
296 |
< |
earr[hp->ns-1] *= 8./3.; |
297 |
< |
earr[(hp->ns-1)*hp->ns] *= 8./3.; |
298 |
< |
earr[(hp->ns-1)*hp->ns + hp->ns-1] *= 8./3.; |
295 |
> |
ep = earr + hp->ns*hp->ns; |
296 |
> |
ep[0] *= 6./3.; |
297 |
> |
ep[hp->ns-1] *= 6./3.; |
298 |
> |
ep[(hp->ns-1)*hp->ns] *= 6./3.; |
299 |
> |
ep[(hp->ns-1)*hp->ns + hp->ns-1] *= 6./3.; |
300 |
|
for (i = 1; i < hp->ns-1; i++) { |
301 |
< |
earr[i*hp->ns] *= 8./5.; |
302 |
< |
earr[i*hp->ns + hp->ns-1] *= 8./5.; |
301 |
> |
ep[i*hp->ns] *= 6./5.; |
302 |
> |
ep[i*hp->ns + hp->ns-1] *= 6./5.; |
303 |
|
} |
304 |
|
for (j = 1; j < hp->ns-1; j++) { |
305 |
< |
earr[j] *= 8./5.; |
306 |
< |
earr[(hp->ns-1)*hp->ns + j] *= 8./5.; |
305 |
> |
ep[j] *= 6./5.; |
306 |
> |
ep[(hp->ns-1)*hp->ns + j] *= 6./5.; |
307 |
|
} |
308 |
+ |
/* blur final map to reduce bias */ |
309 |
+ |
for (i = 0; i < hp->ns-1; i++) { |
310 |
+ |
float *ep2; |
311 |
+ |
ep = earr + i*hp->ns; |
312 |
+ |
ep2 = ep + hp->ns*hp->ns; |
313 |
+ |
for (j = 0; j < hp->ns-1; j++, ep++, ep2++) { |
314 |
+ |
ep[0] += .5*ep2[0] + .125*(ep2[1] + ep2[hp->ns]); |
315 |
+ |
ep[1] += .125*ep2[0]; |
316 |
+ |
ep[hp->ns] += .125*ep2[0]; |
317 |
+ |
} |
318 |
+ |
} |
319 |
|
return(earr); |
320 |
|
} |
321 |
|
|
351 |
|
|
352 |
|
static AMBHEMI * |
353 |
|
samp_hemi( /* sample indirect hemisphere */ |
354 |
< |
COLOR rcol, |
354 |
> |
SCOLOR rcol, |
355 |
|
RAY *r, |
356 |
|
double wt |
357 |
|
) |
358 |
|
{ |
359 |
+ |
int backside = (wt < 0); |
360 |
|
AMBHEMI *hp; |
361 |
|
double d; |
362 |
|
int n, i, j; |
363 |
|
/* insignificance check */ |
364 |
< |
if (bright(rcol) <= FTINY) |
364 |
> |
d = sintens(rcol); |
365 |
> |
if (d <= FTINY) |
366 |
|
return(NULL); |
367 |
|
/* set number of divisions */ |
368 |
+ |
if (backside) wt = -wt; |
369 |
|
if (ambacc <= FTINY && |
370 |
< |
wt > (d = 0.8*intens(rcol)*r->rweight/(ambdiv*minweight))) |
370 |
> |
wt > (d *= 0.8*r->rweight/(ambdiv*minweight + 1e-20))) |
371 |
|
wt = d; /* avoid ray termination */ |
372 |
|
n = sqrt(ambdiv * wt) + 0.5; |
373 |
|
i = 1 + (MINADIV-1)*(ambacc > FTINY); |
377 |
|
hp = (AMBHEMI *)malloc(sizeof(AMBHEMI) + sizeof(AMBSAMP)*(n*n - 1)); |
378 |
|
if (hp == NULL) |
379 |
|
error(SYSTEM, "out of memory in samp_hemi"); |
380 |
+ |
|
381 |
+ |
if (backside) { |
382 |
+ |
hp->atyp = TAMBIENT; |
383 |
+ |
hp->onrm[0] = -r->ron[0]; |
384 |
+ |
hp->onrm[1] = -r->ron[1]; |
385 |
+ |
hp->onrm[2] = -r->ron[2]; |
386 |
+ |
} else { |
387 |
+ |
hp->atyp = RAMBIENT; |
388 |
+ |
VCOPY(hp->onrm, r->ron); |
389 |
+ |
} |
390 |
|
hp->rp = r; |
391 |
|
hp->ns = n; |
392 |
< |
hp->acol[RED] = hp->acol[GRN] = hp->acol[BLU] = 0.0; |
392 |
> |
scolorblack(hp->acol); |
393 |
|
memset(hp->sa, 0, sizeof(AMBSAMP)*n*n); |
394 |
|
hp->sampOK = 0; |
395 |
|
/* assign coefficient */ |
396 |
< |
copycolor(hp->acoef, rcol); |
396 |
> |
copyscolor(hp->acoef, rcol); |
397 |
|
d = 1.0/(n*n); |
398 |
< |
scalecolor(hp->acoef, d); |
398 |
> |
scalescolor(hp->acoef, d); |
399 |
|
/* make tangent plane axes */ |
400 |
< |
if (!getperpendicular(hp->ux, r->ron, 1)) |
400 |
> |
if (!getperpendicular(hp->ux, hp->onrm, 1)) |
401 |
|
error(CONSISTENCY, "bad ray direction in samp_hemi"); |
402 |
< |
VCROSS(hp->uy, r->ron, hp->ux); |
402 |
> |
VCROSS(hp->uy, hp->onrm, hp->ux); |
403 |
|
/* sample divisions */ |
404 |
|
for (i = hp->ns; i--; ) |
405 |
|
for (j = hp->ns; j--; ) |
406 |
|
hp->sampOK += ambsample(hp, i, j, 0); |
407 |
< |
copycolor(rcol, hp->acol); |
407 |
> |
copyscolor(rcol, hp->acol); |
408 |
|
if (!hp->sampOK) { /* utter failure? */ |
409 |
|
free(hp); |
410 |
|
return(NULL); |
416 |
|
if (hp->sampOK <= MINADIV*MINADIV) |
417 |
|
return(hp); /* don't bother super-sampling */ |
418 |
|
n = ambssamp*wt + 0.5; |
419 |
< |
if (n > 8) { /* perform super-sampling? */ |
419 |
> |
if (n >= 4*hp->ns) { /* perform super-sampling? */ |
420 |
|
ambsupersamp(hp, n); |
421 |
< |
copycolor(rcol, hp->acol); |
421 |
> |
copyscolor(rcol, hp->acol); |
422 |
|
} |
423 |
|
return(hp); /* all is well */ |
424 |
|
} |
430 |
|
{ |
431 |
|
if (hp->sa[n1].d <= hp->sa[n2].d) { |
432 |
|
if (hp->sa[n1].d <= hp->sa[n3].d) |
433 |
< |
return(colval(hp->sa[n1].v,CIEY)); |
434 |
< |
return(colval(hp->sa[n3].v,CIEY)); |
433 |
> |
return(hp->sa[n1].v[0]); |
434 |
> |
return(hp->sa[n3].v[0]); |
435 |
|
} |
436 |
|
if (hp->sa[n2].d <= hp->sa[n3].d) |
437 |
< |
return(colval(hp->sa[n2].v,CIEY)); |
438 |
< |
return(colval(hp->sa[n3].v,CIEY)); |
437 |
> |
return(hp->sa[n2].v[0]); |
438 |
> |
return(hp->sa[n3].v[0]); |
439 |
|
} |
440 |
|
|
441 |
|
|
664 |
|
for (j = 0; j < hp->ns-1; j++) { |
665 |
|
comp_fftri(&fftr, hp, AI(hp,0,j), AI(hp,0,j+1)); |
666 |
|
if (hessrow != NULL) |
667 |
< |
comp_hessian(hessrow[j], &fftr, hp->rp->ron); |
667 |
> |
comp_hessian(hessrow[j], &fftr, hp->onrm); |
668 |
|
if (gradrow != NULL) |
669 |
< |
comp_gradient(gradrow[j], &fftr, hp->rp->ron); |
669 |
> |
comp_gradient(gradrow[j], &fftr, hp->onrm); |
670 |
|
} |
671 |
|
/* sum each row of triangles */ |
672 |
|
for (i = 0; i < hp->ns-1; i++) { |
674 |
|
FVECT gradcol; |
675 |
|
comp_fftri(&fftr, hp, AI(hp,i,0), AI(hp,i+1,0)); |
676 |
|
if (hessrow != NULL) |
677 |
< |
comp_hessian(hesscol, &fftr, hp->rp->ron); |
677 |
> |
comp_hessian(hesscol, &fftr, hp->onrm); |
678 |
|
if (gradrow != NULL) |
679 |
< |
comp_gradient(gradcol, &fftr, hp->rp->ron); |
679 |
> |
comp_gradient(gradcol, &fftr, hp->onrm); |
680 |
|
for (j = 0; j < hp->ns-1; j++) { |
681 |
|
FVECT hessdia[3]; /* compute triangle contributions */ |
682 |
|
FVECT graddia; |
686 |
|
/* diagonal (inner) edge */ |
687 |
|
comp_fftri(&fftr, hp, AI(hp,i,j+1), AI(hp,i+1,j)); |
688 |
|
if (hessrow != NULL) { |
689 |
< |
comp_hessian(hessdia, &fftr, hp->rp->ron); |
689 |
> |
comp_hessian(hessdia, &fftr, hp->onrm); |
690 |
|
rev_hessian(hesscol); |
691 |
|
add2hessian(hessian, hessrow[j], hessdia, hesscol, backg); |
692 |
|
} |
693 |
|
if (gradrow != NULL) { |
694 |
< |
comp_gradient(graddia, &fftr, hp->rp->ron); |
694 |
> |
comp_gradient(graddia, &fftr, hp->onrm); |
695 |
|
rev_gradient(gradcol); |
696 |
|
add2gradient(gradient, gradrow[j], graddia, gradcol, backg); |
697 |
|
} |
698 |
|
/* initialize edge in next row */ |
699 |
|
comp_fftri(&fftr, hp, AI(hp,i+1,j+1), AI(hp,i+1,j)); |
700 |
|
if (hessrow != NULL) |
701 |
< |
comp_hessian(hessrow[j], &fftr, hp->rp->ron); |
701 |
> |
comp_hessian(hessrow[j], &fftr, hp->onrm); |
702 |
|
if (gradrow != NULL) |
703 |
< |
comp_gradient(gradrow[j], &fftr, hp->rp->ron); |
703 |
> |
comp_gradient(gradrow[j], &fftr, hp->onrm); |
704 |
|
/* new column edge & paired triangle */ |
705 |
|
backg = back_ambval(hp, AI(hp,i+1,j+1), |
706 |
|
AI(hp,i+1,j), AI(hp,i,j+1)); |
707 |
|
comp_fftri(&fftr, hp, AI(hp,i,j+1), AI(hp,i+1,j+1)); |
708 |
|
if (hessrow != NULL) { |
709 |
< |
comp_hessian(hesscol, &fftr, hp->rp->ron); |
709 |
> |
comp_hessian(hesscol, &fftr, hp->onrm); |
710 |
|
rev_hessian(hessdia); |
711 |
|
add2hessian(hessian, hessrow[j], hessdia, hesscol, backg); |
712 |
|
if (i < hp->ns-2) |
713 |
|
rev_hessian(hessrow[j]); |
714 |
|
} |
715 |
|
if (gradrow != NULL) { |
716 |
< |
comp_gradient(gradcol, &fftr, hp->rp->ron); |
716 |
> |
comp_gradient(gradcol, &fftr, hp->onrm); |
717 |
|
rev_gradient(graddia); |
718 |
|
add2gradient(gradient, gradrow[j], graddia, gradcol, backg); |
719 |
|
if (i < hp->ns-2) |
749 |
|
/* use vector for azimuth + 90deg */ |
750 |
|
VSUB(vd, ap->p, hp->rp->rop); |
751 |
|
/* brightness over cosine factor */ |
752 |
< |
gfact = colval(ap->v,CIEY) / DOT(hp->rp->ron, vd); |
752 |
> |
gfact = ap->v[0] / DOT(hp->onrm, vd); |
753 |
|
/* sine = proj_radius/vd_length */ |
754 |
|
dgsum[0] -= DOT(uv[1], vd) * gfact; |
755 |
|
dgsum[1] += DOT(uv[0], vd) * gfact; |
805 |
|
|
806 |
|
int |
807 |
|
doambient( /* compute ambient component */ |
808 |
< |
COLOR rcol, /* input/output color */ |
808 |
> |
SCOLOR rcol, /* input/output color */ |
809 |
|
RAY *r, |
810 |
< |
double wt, |
810 |
> |
double wt, /* negative for back side */ |
811 |
|
FVECT uv[2], /* returned (optional) */ |
812 |
|
float ra[2], /* returned (optional) */ |
813 |
|
float pg[2], /* returned (optional) */ |
839 |
|
free(hp); /* Hessian not requested/possible */ |
840 |
|
return(-1); /* value-only return value */ |
841 |
|
} |
842 |
< |
if ((d = bright(rcol)) > FTINY) { /* normalize Y values */ |
843 |
< |
d = 0.99*(hp->ns*hp->ns)/d; |
842 |
> |
if ((d = scolor_mean(rcol)) > FTINY) { |
843 |
> |
d = 0.99*(hp->ns*hp->ns)/d; /* normalize avg. values */ |
844 |
|
K = 0.01; |
845 |
|
} else { /* or fall back on geometric Hessian */ |
846 |
|
K = 1.0; |
848 |
|
dg = NULL; |
849 |
|
crlp = NULL; |
850 |
|
} |
851 |
< |
ap = hp->sa; /* relative Y channel from here on... */ |
851 |
> |
ap = hp->sa; /* single channel from here on... */ |
852 |
|
for (i = hp->ns*hp->ns; i--; ap++) |
853 |
< |
colval(ap->v,CIEY) = bright(ap->v)*d + K; |
853 |
> |
ap->v[0] = scolor_mean(ap->v)*d + K; |
854 |
|
|
855 |
|
if (uv == NULL) /* make sure we have axis pointers */ |
856 |
|
uv = my_uv; |
874 |
|
if (ra[1] < minarad) |
875 |
|
ra[1] = minarad; |
876 |
|
} |
877 |
< |
ra[0] *= d = 1.0/sqrt(wt); |
877 |
> |
ra[0] *= d = 1.0/sqrt(fabs(wt)); |
878 |
|
if ((ra[1] *= d) > 2.0*ra[0]) |
879 |
|
ra[1] = 2.0*ra[0]; |
880 |
|
if (ra[1] > maxarad) { |
883 |
|
ra[0] = maxarad; |
884 |
|
} |
885 |
|
/* flag encroached directions */ |
886 |
< |
if (crlp != NULL) |
886 |
> |
if (crlp != NULL) /* XXX doesn't update with changes to ambacc */ |
887 |
|
*crlp = ambcorral(hp, uv, ra[0]*ambacc, ra[1]*ambacc); |
888 |
|
if (pg != NULL) { /* cap gradient if necessary */ |
889 |
|
d = pg[0]*pg[0]*ra[0]*ra[0] + pg[1]*pg[1]*ra[1]*ra[1]; |
897 |
|
free(hp); /* clean up and return */ |
898 |
|
return(1); |
899 |
|
} |
782 |
– |
|