| 11 | 
  | 
#include  "face.h" | 
| 12 | 
  | 
#include  "cone.h" | 
| 13 | 
  | 
#include  "source.h" | 
| 14 | 
+ | 
#include  "paths.h" | 
| 15 | 
  | 
 | 
| 16 | 
< | 
#ifndef NBSDFSAMPS | 
| 17 | 
< | 
#define NBSDFSAMPS      256             /* BSDF resampling count */ | 
| 16 | 
> | 
#ifndef R_EPS | 
| 17 | 
> | 
#define R_EPS           0.005           /* relative epsilon for ray origin */ | 
| 18 | 
  | 
#endif | 
| 19 | 
  | 
 | 
| 20 | 
  | 
COLORV *        distarr = NULL;         /* distribution array */ | 
| 21 | 
  | 
int             distsiz = 0; | 
| 21 | 
– | 
COLORV *        direct_discount = NULL; /* amount to take off direct */ | 
| 22 | 
  | 
 | 
| 23 | 
  | 
 | 
| 24 | 
  | 
void | 
| 28 | 
  | 
{ | 
| 29 | 
  | 
        if (siz <= 0) { | 
| 30 | 
  | 
                if (distsiz > 0) | 
| 31 | 
< | 
                        free((void *)distarr); | 
| 31 | 
> | 
                        free(distarr); | 
| 32 | 
  | 
                distarr = NULL; | 
| 33 | 
  | 
                distsiz = 0; | 
| 34 | 
  | 
                return; | 
| 35 | 
  | 
        } | 
| 36 | 
  | 
        if (distsiz < siz) { | 
| 37 | 
  | 
                if (distsiz > 0) | 
| 38 | 
< | 
                        free((void *)distarr); | 
| 38 | 
> | 
                        free(distarr); | 
| 39 | 
  | 
                distarr = (COLORV *)malloc(sizeof(COLOR)*siz); | 
| 40 | 
  | 
                if (distarr == NULL) | 
| 41 | 
  | 
                        error(SYSTEM, "out of memory in newdist"); | 
| 45 | 
  | 
} | 
| 46 | 
  | 
 | 
| 47 | 
  | 
 | 
| 48 | 
– | 
static void | 
| 49 | 
– | 
new_discount()                  /* allocate space for direct contrib. record */ | 
| 50 | 
– | 
{ | 
| 51 | 
– | 
        if (distsiz <= 0) | 
| 52 | 
– | 
                return; | 
| 53 | 
– | 
        direct_discount = (COLORV *)calloc(distsiz, sizeof(COLOR)); | 
| 54 | 
– | 
        if (direct_discount == NULL) | 
| 55 | 
– | 
                error(SYSTEM, "out of memory in new_discount"); | 
| 56 | 
– | 
} | 
| 57 | 
– | 
 | 
| 58 | 
– | 
 | 
| 59 | 
– | 
static void | 
| 60 | 
– | 
done_discount()                 /* clear off direct contrib. record */ | 
| 61 | 
– | 
{ | 
| 62 | 
– | 
        if (direct_discount == NULL) | 
| 63 | 
– | 
                return; | 
| 64 | 
– | 
        free((void *)direct_discount); | 
| 65 | 
– | 
        direct_discount = NULL; | 
| 66 | 
– | 
} | 
| 67 | 
– | 
 | 
| 68 | 
– | 
 | 
| 48 | 
  | 
int | 
| 49 | 
  | 
process_ray(                    /* process a ray result or report error */ | 
| 50 | 
  | 
        RAY *r, | 
| 62 | 
  | 
        multcolor(r->rcol, r->rcoef);   /* in case it's a source ray */ | 
| 63 | 
  | 
        colp = &distarr[r->rno * 3]; | 
| 64 | 
  | 
        addcolor(colp, r->rcol); | 
| 86 | 
– | 
        if (r->rsrc >= 0 &&             /* remember source contrib. */ | 
| 87 | 
– | 
                        direct_discount != NULL) { | 
| 88 | 
– | 
                colp = &direct_discount[r->rno * 3]; | 
| 89 | 
– | 
                addcolor(colp, r->rcol); | 
| 90 | 
– | 
        } | 
| 65 | 
  | 
        return(1); | 
| 66 | 
  | 
} | 
| 67 | 
  | 
 | 
| 81 | 
  | 
        VCOPY(myRay.rorg, org); | 
| 82 | 
  | 
        VCOPY(myRay.rdir, dir); | 
| 83 | 
  | 
        myRay.rmax = .0; | 
| 84 | 
< | 
        rayorigin(&myRay, PRIMARY, NULL, NULL); | 
| 84 | 
> | 
        rayorigin(&myRay, PRIMARY|SPECULAR, NULL, NULL); | 
| 85 | 
  | 
        myRay.rno = ndx; | 
| 86 | 
  | 
                                        /* queue ray, check result */ | 
| 87 | 
  | 
        process_ray(&myRay, ray_pqueue(&myRay)); | 
| 92 | 
  | 
srcsamps(                       /* sample sources from this surface position */ | 
| 93 | 
  | 
        struct illum_args *il, | 
| 94 | 
  | 
        FVECT org, | 
| 95 | 
< | 
        FVECT nrm, | 
| 95 | 
> | 
        double eps, | 
| 96 | 
  | 
        MAT4 ixfm | 
| 97 | 
  | 
) | 
| 98 | 
  | 
{ | 
| 99 | 
< | 
        int  nalt, nazi; | 
| 99 | 
> | 
        int  nalt=1, nazi=1; | 
| 100 | 
  | 
        SRCINDEX  si; | 
| 101 | 
  | 
        RAY  sr; | 
| 102 | 
  | 
        FVECT   v; | 
| 103 | 
  | 
        double  d; | 
| 104 | 
  | 
        int  i, j; | 
| 105 | 
  | 
                                                /* get sampling density */ | 
| 106 | 
< | 
        if (il->sampdens <= 0) { | 
| 133 | 
< | 
                nalt = nazi = 1; | 
| 134 | 
< | 
        } else { | 
| 106 | 
> | 
        if (il->sampdens > 0) { | 
| 107 | 
  | 
                i = PI * il->sampdens; | 
| 108 | 
  | 
                nalt = sqrt(i/PI) + .5; | 
| 109 | 
  | 
                nazi = PI*nalt + .5; | 
| 111 | 
  | 
        initsrcindex(&si);                      /* loop over (sub)sources */ | 
| 112 | 
  | 
        for ( ; ; ) { | 
| 113 | 
  | 
                VCOPY(sr.rorg, org);            /* pick side to shoot from */ | 
| 142 | 
– | 
                if (il->sd != NULL) { | 
| 143 | 
– | 
                        int  sn = si.sn; | 
| 144 | 
– | 
                        if (si.sp+1 >= si.np) ++sn; | 
| 145 | 
– | 
                        if (sn >= nsources) break; | 
| 146 | 
– | 
                        if (source[sn].sflags & SDISTANT) | 
| 147 | 
– | 
                                d = DOT(source[sn].sloc, nrm); | 
| 148 | 
– | 
                        else { | 
| 149 | 
– | 
                                VSUB(v, source[sn].sloc, org); | 
| 150 | 
– | 
                                d = DOT(v, nrm); | 
| 151 | 
– | 
                        } | 
| 152 | 
– | 
                } else | 
| 153 | 
– | 
                        d = 1.0;                /* only transmission */ | 
| 154 | 
– | 
                if (d < 0.0) | 
| 155 | 
– | 
                        d = -1.0001*il->thick - 5.*FTINY; | 
| 156 | 
– | 
                else | 
| 157 | 
– | 
                        d = 5.*FTINY; | 
| 158 | 
– | 
                for (i = 3; i--; ) | 
| 159 | 
– | 
                        sr.rorg[i] += d*nrm[i]; | 
| 114 | 
  | 
                samplendx++;                    /* increment sample counter */ | 
| 115 | 
  | 
                if (!srcray(&sr, NULL, &si)) | 
| 116 | 
  | 
                        break;                  /* end of sources */ | 
| 119 | 
  | 
                        multv3(v, sr.rdir, ixfm); | 
| 120 | 
  | 
                else | 
| 121 | 
  | 
                        VCOPY(v, sr.rdir); | 
| 122 | 
< | 
                if (il->sd != NULL) { | 
| 123 | 
< | 
                        i = getBSDF_incndx(il->sd, v); | 
| 124 | 
< | 
                        if (i < 0) | 
| 125 | 
< | 
                                continue;       /* must not be important */ | 
| 126 | 
< | 
                        sr.rno = i; | 
| 173 | 
< | 
                        d = 1.0/getBSDF_incohm(il->sd, i); | 
| 174 | 
< | 
                } else { | 
| 175 | 
< | 
                        if (v[2] >= -FTINY) | 
| 176 | 
< | 
                                continue;       /* only sample transmission */ | 
| 177 | 
< | 
                        v[0] = -v[0]; v[1] = -v[1]; v[2] = -v[2]; | 
| 178 | 
< | 
                        sr.rno = flatindex(v, nalt, nazi); | 
| 179 | 
< | 
                        d = nalt*nazi*(1./PI) * v[2]; | 
| 180 | 
< | 
                } | 
| 122 | 
> | 
                if (v[2] >= -FTINY) | 
| 123 | 
> | 
                        continue;               /* only sample transmission */ | 
| 124 | 
> | 
                v[0] = -v[0]; v[1] = -v[1]; v[2] = -v[2]; | 
| 125 | 
> | 
                sr.rno = flatindex(v, nalt, nazi); | 
| 126 | 
> | 
                d = nalt*nazi*(1./PI) * v[2]; | 
| 127 | 
  | 
                d *= si.dom;                    /* solid angle correction */ | 
| 128 | 
  | 
                scalecolor(sr.rcoef, d); | 
| 129 | 
+ | 
                VSUM(sr.rorg, sr.rorg, sr.rdir, -eps); | 
| 130 | 
  | 
                process_ray(&sr, ray_pqueue(&sr)); | 
| 131 | 
  | 
        } | 
| 132 | 
  | 
} | 
| 149 | 
  | 
        FVECT  n | 
| 150 | 
  | 
) | 
| 151 | 
  | 
{ | 
| 152 | 
< | 
        register int  i; | 
| 206 | 
< | 
 | 
| 207 | 
< | 
        v[0] = v[1] = v[2] = 0.0; | 
| 208 | 
< | 
        for (i = 0; i < 3; i++) | 
| 209 | 
< | 
                if (n[i] < 0.6 && n[i] > -0.6) | 
| 210 | 
< | 
                        break; | 
| 211 | 
< | 
        v[i] = 1.0; | 
| 212 | 
< | 
        fcross(u, v, n); | 
| 213 | 
< | 
        normalize(u); | 
| 152 | 
> | 
        getperpendicular(u, n, 1); | 
| 153 | 
  | 
        fcross(v, n, u); | 
| 154 | 
  | 
} | 
| 155 | 
  | 
 | 
| 156 | 
  | 
 | 
| 157 | 
  | 
static void | 
| 158 | 
  | 
rounddir(               /* compute uniform spherical direction */ | 
| 159 | 
< | 
        register FVECT  dv, | 
| 159 | 
> | 
        FVECT  dv, | 
| 160 | 
  | 
        double  alt, | 
| 161 | 
  | 
        double  azi | 
| 162 | 
  | 
) | 
| 187 | 
  | 
        dv[2] = sqrt(1. - alt); | 
| 188 | 
  | 
} | 
| 189 | 
  | 
 | 
| 190 | 
+ | 
 | 
| 191 | 
  | 
int | 
| 192 | 
  | 
flatindex(              /* compute index for hemispherical direction */ | 
| 193 | 
  | 
        FVECT   dv, | 
| 237 | 
  | 
        FVECT  dn, org, dir; | 
| 238 | 
  | 
        FVECT  u, v; | 
| 239 | 
  | 
        double  ur[2], vr[2]; | 
| 240 | 
+ | 
        double  epsilon; | 
| 241 | 
  | 
        MAT4  xfm; | 
| 242 | 
+ | 
        char  xfrot[64]; | 
| 243 | 
  | 
        int  nallow; | 
| 244 | 
  | 
        FACE  *fa; | 
| 245 | 
  | 
        int  i, j; | 
| 250 | 
  | 
                return(my_default(ob, il, nm)); | 
| 251 | 
  | 
        } | 
| 252 | 
  | 
                                /* set up sampling */ | 
| 253 | 
< | 
        if (il->sd != NULL) { | 
| 254 | 
< | 
                if (!getBSDF_xfm(xfm, fa->norm, il->udir)) { | 
| 313 | 
< | 
                        objerror(ob, WARNING, "illegal up direction"); | 
| 314 | 
< | 
                        freeface(ob); | 
| 315 | 
< | 
                        return(my_default(ob, il, nm)); | 
| 316 | 
< | 
                } | 
| 317 | 
< | 
                n = il->sd->ninc; | 
| 253 | 
> | 
        if (il->sampdens <= 0) { | 
| 254 | 
> | 
                nalt = nazi = 1;        /* diffuse assumption */ | 
| 255 | 
  | 
        } else { | 
| 256 | 
< | 
                if (il->sampdens <= 0) { | 
| 257 | 
< | 
                        nalt = nazi = 1;        /* diffuse assumption */ | 
| 258 | 
< | 
                } else { | 
| 322 | 
< | 
                        n = PI * il->sampdens; | 
| 323 | 
< | 
                        nalt = sqrt(n/PI) + .5; | 
| 324 | 
< | 
                        nazi = PI*nalt + .5; | 
| 325 | 
< | 
                } | 
| 326 | 
< | 
                n = nazi*nalt; | 
| 256 | 
> | 
                n = PI * il->sampdens; | 
| 257 | 
> | 
                nalt = sqrt(n/PI) + .5; | 
| 258 | 
> | 
                nazi = PI*nalt + .5; | 
| 259 | 
  | 
        } | 
| 260 | 
+ | 
        n = nazi*nalt; | 
| 261 | 
  | 
        newdist(n); | 
| 262 | 
  | 
                                /* take first edge >= sqrt(area) */ | 
| 263 | 
  | 
        for (j = fa->nv-1, i = 0; i < fa->nv; j = i++) { | 
| 287 | 
  | 
        dim[0] = random(); | 
| 288 | 
  | 
                                /* sample polygon */ | 
| 289 | 
  | 
        nallow = 5*n*il->nsamps; | 
| 290 | 
+ | 
        epsilon = R_EPS*sqrt(fa->area); | 
| 291 | 
  | 
        for (dim[1] = 0; dim[1] < n; dim[1]++) | 
| 292 | 
  | 
                for (i = 0; i < il->nsamps; i++) { | 
| 293 | 
  | 
                                        /* randomize direction */ | 
| 294 | 
  | 
                    h = ilhash(dim, 2) + i; | 
| 295 | 
< | 
                    if (il->sd != NULL) { | 
| 296 | 
< | 
                        r_BSDF_incvec(dir, il->sd, dim[1], urand(h), xfm); | 
| 297 | 
< | 
                    } else { | 
| 298 | 
< | 
                        multisamp(sp, 2, urand(h)); | 
| 299 | 
< | 
                        alti = dim[1]/nazi; | 
| 300 | 
< | 
                        r1 = (alti + sp[0])/nalt; | 
| 367 | 
< | 
                        r2 = (dim[1] - alti*nazi + sp[1] - .5)/nazi; | 
| 368 | 
< | 
                        flatdir(dn, r1, r2); | 
| 369 | 
< | 
                        for (j = 0; j < 3; j++) | 
| 295 | 
> | 
                    multisamp(sp, 2, urand(h)); | 
| 296 | 
> | 
                    alti = dim[1]/nazi; | 
| 297 | 
> | 
                    r1 = (alti + sp[0])/nalt; | 
| 298 | 
> | 
                    r2 = (dim[1] - alti*nazi + sp[1] - .5)/nazi; | 
| 299 | 
> | 
                    flatdir(dn, r1, r2); | 
| 300 | 
> | 
                    for (j = 0; j < 3; j++) | 
| 301 | 
  | 
                            dir[j] = -dn[0]*u[j] - dn[1]*v[j] - | 
| 302 | 
  | 
                                                dn[2]*fa->norm[j]; | 
| 372 | 
– | 
                    } | 
| 303 | 
  | 
                                        /* randomize location */ | 
| 304 | 
  | 
                    do { | 
| 305 | 
  | 
                        multisamp(sp, 2, urand(h+4862+nallow)); | 
| 315 | 
  | 
                        freeface(ob); | 
| 316 | 
  | 
                        return(my_default(ob, il, nm)); | 
| 317 | 
  | 
                    } | 
| 318 | 
< | 
                    if (il->sd != NULL && DOT(dir, fa->norm) < -FTINY) | 
| 389 | 
< | 
                        r1 = -1.0001*il->thick - 5.*FTINY; | 
| 390 | 
< | 
                    else | 
| 391 | 
< | 
                        r1 = 5.*FTINY; | 
| 392 | 
< | 
                    for (j = 0; j < 3; j++) | 
| 393 | 
< | 
                        org[j] += r1*fa->norm[j]; | 
| 318 | 
> | 
                    VSUM(org, org, dir, -epsilon); | 
| 319 | 
  | 
                                        /* send sample */ | 
| 320 | 
  | 
                    raysamp(dim[1], org, dir); | 
| 321 | 
  | 
                } | 
| 322 | 
  | 
                                /* add in direct component? */ | 
| 323 | 
< | 
        if (!directvis && (il->flags & IL_LIGHT || il->sd != NULL)) { | 
| 323 | 
> | 
        if (il->flags & IL_LIGHT) { | 
| 324 | 
  | 
                MAT4    ixfm; | 
| 325 | 
< | 
                if (il->sd == NULL) { | 
| 326 | 
< | 
                        for (i = 3; i--; ) { | 
| 327 | 
< | 
                                ixfm[i][0] = u[i]; | 
| 328 | 
< | 
                                ixfm[i][1] = v[i]; | 
| 329 | 
< | 
                                ixfm[i][2] = fa->norm[i]; | 
| 405 | 
< | 
                                ixfm[i][3] = 0.; | 
| 406 | 
< | 
                        } | 
| 407 | 
< | 
                        ixfm[3][0] = ixfm[3][1] = ixfm[3][2] = 0.; | 
| 408 | 
< | 
                        ixfm[3][3] = 1.; | 
| 409 | 
< | 
                } else { | 
| 410 | 
< | 
                        if (!invmat4(ixfm, xfm)) | 
| 411 | 
< | 
                                objerror(ob, INTERNAL, | 
| 412 | 
< | 
                                        "cannot invert BSDF transform"); | 
| 413 | 
< | 
                        if (!(il->flags & IL_LIGHT)) | 
| 414 | 
< | 
                                new_discount(); | 
| 325 | 
> | 
                for (i = 3; i--; ) { | 
| 326 | 
> | 
                        ixfm[i][0] = u[i]; | 
| 327 | 
> | 
                        ixfm[i][1] = v[i]; | 
| 328 | 
> | 
                        ixfm[i][2] = fa->norm[i]; | 
| 329 | 
> | 
                        ixfm[i][3] = 0.; | 
| 330 | 
  | 
                } | 
| 331 | 
+ | 
                ixfm[3][0] = ixfm[3][1] = ixfm[3][2] = 0.; | 
| 332 | 
+ | 
                ixfm[3][3] = 1.; | 
| 333 | 
  | 
                dim[0] = random(); | 
| 334 | 
  | 
                nallow = 10*il->nsamps; | 
| 335 | 
  | 
                for (i = 0; i < il->nsamps; i++) { | 
| 350 | 
  | 
                        return(my_default(ob, il, nm)); | 
| 351 | 
  | 
                    } | 
| 352 | 
  | 
                                        /* sample source rays */ | 
| 353 | 
< | 
                    srcsamps(il, org, fa->norm, ixfm); | 
| 353 | 
> | 
                    srcsamps(il, org, epsilon, ixfm); | 
| 354 | 
  | 
                } | 
| 355 | 
  | 
        } | 
| 356 | 
  | 
                                /* wait for all rays to finish */ | 
| 357 | 
  | 
        rayclean(); | 
| 441 | 
– | 
        if (il->sd != NULL) {   /* run distribution through BSDF */ | 
| 442 | 
– | 
                nalt = sqrt(il->sd->nout/PI) + .5; | 
| 443 | 
– | 
                nazi = PI*nalt + .5; | 
| 444 | 
– | 
                redistribute(il->sd, nalt, nazi, u, v, fa->norm, xfm); | 
| 445 | 
– | 
                done_discount(); | 
| 446 | 
– | 
                if (!il->sampdens) | 
| 447 | 
– | 
                        il->sampdens = nalt*nazi/PI + .999; | 
| 448 | 
– | 
        } | 
| 358 | 
  | 
                                /* write out the face and its distribution */ | 
| 359 | 
  | 
        if (average(il, distarr, n)) { | 
| 360 | 
  | 
                if (il->sampdens > 0) | 
| 370 | 
  | 
 | 
| 371 | 
  | 
int | 
| 372 | 
  | 
my_sphere(      /* make an illum sphere */ | 
| 373 | 
< | 
        register OBJREC  *ob, | 
| 373 | 
> | 
        OBJREC  *ob, | 
| 374 | 
  | 
        struct illum_args  *il, | 
| 375 | 
  | 
        char  *nm | 
| 376 | 
  | 
) | 
| 380 | 
  | 
        double  sp[4], r1, r2, r3; | 
| 381 | 
  | 
        FVECT  org, dir; | 
| 382 | 
  | 
        FVECT  u, v; | 
| 383 | 
< | 
        register int  i, j; | 
| 383 | 
> | 
        int  i, j; | 
| 384 | 
  | 
                                /* check arguments */ | 
| 385 | 
  | 
        if (ob->oargs.nfargs != 4) | 
| 386 | 
  | 
                objerror(ob, USER, "bad # of arguments"); | 
| 392 | 
  | 
                nalt = sqrt(2./PI*n) + .5; | 
| 393 | 
  | 
                nazi = PI/2.*nalt + .5; | 
| 394 | 
  | 
        } | 
| 486 | 
– | 
        if (il->sd != NULL) | 
| 487 | 
– | 
                objerror(ob, WARNING, "BSDF ignored"); | 
| 395 | 
  | 
        n = nalt*nazi; | 
| 396 | 
  | 
        newdist(n); | 
| 397 | 
  | 
        dim[0] = random(); | 
| 446 | 
  | 
        int  dim[2]; | 
| 447 | 
  | 
        int  n, nalt, nazi, alti; | 
| 448 | 
  | 
        double  sp[2], r1, r2, r3; | 
| 449 | 
+ | 
        double  epsilon; | 
| 450 | 
  | 
        int  h; | 
| 451 | 
  | 
        FVECT  dn, org, dir; | 
| 452 | 
  | 
        FVECT  u, v; | 
| 453 | 
  | 
        MAT4  xfm; | 
| 454 | 
  | 
        CONE  *co; | 
| 455 | 
  | 
        int  i, j; | 
| 456 | 
< | 
                                /* get/check arguments */ | 
| 456 | 
> | 
                                        /* get/check arguments */ | 
| 457 | 
  | 
        co = getcone(ob, 0); | 
| 458 | 
< | 
                                /* set up sampling */ | 
| 459 | 
< | 
        if (il->sd != NULL) { | 
| 460 | 
< | 
                if (!getBSDF_xfm(xfm, co->ad, il->udir)) { | 
| 461 | 
< | 
                        objerror(ob, WARNING, "illegal up direction"); | 
| 462 | 
< | 
                        freecone(ob); | 
| 555 | 
< | 
                        return(my_default(ob, il, nm)); | 
| 556 | 
< | 
                } | 
| 557 | 
< | 
                n = il->sd->ninc; | 
| 458 | 
> | 
        if (co == NULL) | 
| 459 | 
> | 
                objerror(ob, USER, "cannot create illum"); | 
| 460 | 
> | 
                                        /* set up sampling */ | 
| 461 | 
> | 
        if (il->sampdens <= 0) { | 
| 462 | 
> | 
                nalt = nazi = 1;        /* diffuse assumption */ | 
| 463 | 
  | 
        } else { | 
| 464 | 
< | 
                if (il->sampdens <= 0) { | 
| 465 | 
< | 
                        nalt = nazi = 1;        /* diffuse assumption */ | 
| 466 | 
< | 
                } else { | 
| 562 | 
< | 
                        n = PI * il->sampdens; | 
| 563 | 
< | 
                        nalt = sqrt(n/PI) + .5; | 
| 564 | 
< | 
                        nazi = PI*nalt + .5; | 
| 565 | 
< | 
                } | 
| 566 | 
< | 
                n = nazi*nalt; | 
| 464 | 
> | 
                n = PI * il->sampdens; | 
| 465 | 
> | 
                nalt = sqrt(n/PI) + .5; | 
| 466 | 
> | 
                nazi = PI*nalt + .5; | 
| 467 | 
  | 
        } | 
| 468 | 
+ | 
        epsilon = R_EPS*CO_R1(co); | 
| 469 | 
+ | 
        n = nazi*nalt; | 
| 470 | 
  | 
        newdist(n); | 
| 471 | 
  | 
        mkaxes(u, v, co->ad); | 
| 472 | 
  | 
        dim[0] = random(); | 
| 476 | 
  | 
                                        /* next sample point */ | 
| 477 | 
  | 
                    h = ilhash(dim,2) + i; | 
| 478 | 
  | 
                                        /* randomize direction */ | 
| 479 | 
< | 
                    if (il->sd != NULL) { | 
| 480 | 
< | 
                        r_BSDF_incvec(dir, il->sd, dim[1], urand(h), xfm); | 
| 481 | 
< | 
                    } else { | 
| 482 | 
< | 
                        multisamp(sp, 2, urand(h)); | 
| 483 | 
< | 
                        alti = dim[1]/nazi; | 
| 484 | 
< | 
                        r1 = (alti + sp[0])/nalt; | 
| 485 | 
< | 
                        r2 = (dim[1] - alti*nazi + sp[1] - .5)/nazi; | 
| 584 | 
< | 
                        flatdir(dn, r1, r2); | 
| 585 | 
< | 
                        for (j = 0; j < 3; j++) | 
| 586 | 
< | 
                                dir[j] = -dn[0]*u[j] - dn[1]*v[j] - dn[2]*co->ad[j]; | 
| 587 | 
< | 
                    } | 
| 479 | 
> | 
                    multisamp(sp, 2, urand(h)); | 
| 480 | 
> | 
                    alti = dim[1]/nazi; | 
| 481 | 
> | 
                    r1 = (alti + sp[0])/nalt; | 
| 482 | 
> | 
                    r2 = (dim[1] - alti*nazi + sp[1] - .5)/nazi; | 
| 483 | 
> | 
                    flatdir(dn, r1, r2); | 
| 484 | 
> | 
                    for (j = 0; j < 3; j++) | 
| 485 | 
> | 
                        dir[j] = -dn[0]*u[j] - dn[1]*v[j] - dn[2]*co->ad[j]; | 
| 486 | 
  | 
                                        /* randomize location */ | 
| 487 | 
  | 
                    multisamp(sp, 2, urand(h+8371)); | 
| 488 | 
  | 
                    r3 = sqrt(CO_R0(co)*CO_R0(co) + | 
| 490 | 
  | 
                    r2 = 2.*PI*sp[1]; | 
| 491 | 
  | 
                    r1 = r3*cos(r2); | 
| 492 | 
  | 
                    r2 = r3*sin(r2); | 
| 595 | 
– | 
                    if (il->sd != NULL && DOT(dir, co->ad) < -FTINY) | 
| 596 | 
– | 
                        r3 = -1.0001*il->thick - 5.*FTINY; | 
| 597 | 
– | 
                    else | 
| 598 | 
– | 
                        r3 = 5.*FTINY; | 
| 493 | 
  | 
                    for (j = 0; j < 3; j++) | 
| 494 | 
  | 
                        org[j] = CO_P0(co)[j] + r1*u[j] + r2*v[j] + | 
| 495 | 
< | 
                                                r3*co->ad[j]; | 
| 495 | 
> | 
                                                epsilon*co->ad[j]; | 
| 496 | 
  | 
                                        /* send sample */ | 
| 497 | 
  | 
                    raysamp(dim[1], org, dir); | 
| 498 | 
  | 
                } | 
| 499 | 
  | 
                                /* add in direct component? */ | 
| 500 | 
< | 
        if (!directvis && (il->flags & IL_LIGHT || il->sd != NULL)) { | 
| 500 | 
> | 
        if (il->flags & IL_LIGHT) { | 
| 501 | 
  | 
                MAT4    ixfm; | 
| 502 | 
< | 
                if (il->sd == NULL) { | 
| 503 | 
< | 
                        for (i = 3; i--; ) { | 
| 504 | 
< | 
                                ixfm[i][0] = u[i]; | 
| 505 | 
< | 
                                ixfm[i][1] = v[i]; | 
| 506 | 
< | 
                                ixfm[i][2] = co->ad[i]; | 
| 613 | 
< | 
                                ixfm[i][3] = 0.; | 
| 614 | 
< | 
                        } | 
| 615 | 
< | 
                        ixfm[3][0] = ixfm[3][1] = ixfm[3][2] = 0.; | 
| 616 | 
< | 
                        ixfm[3][3] = 1.; | 
| 617 | 
< | 
                } else { | 
| 618 | 
< | 
                        if (!invmat4(ixfm, xfm)) | 
| 619 | 
< | 
                                objerror(ob, INTERNAL, | 
| 620 | 
< | 
                                        "cannot invert BSDF transform"); | 
| 621 | 
< | 
                        if (!(il->flags & IL_LIGHT)) | 
| 622 | 
< | 
                                new_discount(); | 
| 502 | 
> | 
                for (i = 3; i--; ) { | 
| 503 | 
> | 
                        ixfm[i][0] = u[i]; | 
| 504 | 
> | 
                        ixfm[i][1] = v[i]; | 
| 505 | 
> | 
                        ixfm[i][2] = co->ad[i]; | 
| 506 | 
> | 
                        ixfm[i][3] = 0.; | 
| 507 | 
  | 
                } | 
| 508 | 
+ | 
                ixfm[3][0] = ixfm[3][1] = ixfm[3][2] = 0.; | 
| 509 | 
+ | 
                ixfm[3][3] = 1.; | 
| 510 | 
  | 
                dim[0] = random(); | 
| 511 | 
  | 
                for (i = 0; i < il->nsamps; i++) { | 
| 512 | 
  | 
                                        /* randomize location */ | 
| 520 | 
  | 
                    for (j = 0; j < 3; j++) | 
| 521 | 
  | 
                        org[j] = CO_P0(co)[j] + r1*u[j] + r2*v[j]; | 
| 522 | 
  | 
                                        /* sample source rays */ | 
| 523 | 
< | 
                    srcsamps(il, org, co->ad, ixfm); | 
| 523 | 
> | 
                    srcsamps(il, org, epsilon, ixfm); | 
| 524 | 
  | 
                } | 
| 525 | 
  | 
        } | 
| 526 | 
  | 
                                /* wait for all rays to finish */ | 
| 527 | 
  | 
        rayclean(); | 
| 642 | 
– | 
        if (il->sd != NULL) {   /* run distribution through BSDF */ | 
| 643 | 
– | 
                nalt = sqrt(il->sd->nout/PI) + .5; | 
| 644 | 
– | 
                nazi = PI*nalt + .5; | 
| 645 | 
– | 
                redistribute(il->sd, nalt, nazi, u, v, co->ad, xfm); | 
| 646 | 
– | 
                done_discount(); | 
| 647 | 
– | 
                if (!il->sampdens) | 
| 648 | 
– | 
                        il->sampdens = nalt*nazi/PI + .999; | 
| 649 | 
– | 
        } | 
| 528 | 
  | 
                                /* write out the ring and its distribution */ | 
| 529 | 
  | 
        if (average(il, distarr, n)) { | 
| 530 | 
  | 
                if (il->sampdens > 0) | 
| 535 | 
  | 
                                /* clean up */ | 
| 536 | 
  | 
        freecone(ob); | 
| 537 | 
  | 
        return(1); | 
| 660 | 
– | 
} | 
| 661 | 
– | 
 | 
| 662 | 
– | 
 | 
| 663 | 
– | 
void | 
| 664 | 
– | 
redistribute(           /* pass distarr ray sums through BSDF */ | 
| 665 | 
– | 
        struct BSDF_data *b, | 
| 666 | 
– | 
        int nalt, | 
| 667 | 
– | 
        int nazi, | 
| 668 | 
– | 
        FVECT u, | 
| 669 | 
– | 
        FVECT v, | 
| 670 | 
– | 
        FVECT w, | 
| 671 | 
– | 
        MAT4 xm | 
| 672 | 
– | 
) | 
| 673 | 
– | 
{ | 
| 674 | 
– | 
        int     nout = 0; | 
| 675 | 
– | 
        MAT4    mymat, inmat; | 
| 676 | 
– | 
        COLORV  *idist; | 
| 677 | 
– | 
        COLORV  *cp; | 
| 678 | 
– | 
        FVECT   dv; | 
| 679 | 
– | 
        double  wt; | 
| 680 | 
– | 
        int     i, j, k, c, o; | 
| 681 | 
– | 
        COLOR   col, cinc; | 
| 682 | 
– | 
                                        /* copy incoming distribution */ | 
| 683 | 
– | 
        if (b->ninc > distsiz) | 
| 684 | 
– | 
                error(INTERNAL, "error 1 in redistribute"); | 
| 685 | 
– | 
        idist = (COLORV *)malloc(sizeof(COLOR)*b->ninc); | 
| 686 | 
– | 
        if (idist == NULL) | 
| 687 | 
– | 
                error(SYSTEM, "out of memory in redistribute"); | 
| 688 | 
– | 
        memcpy(idist, distarr, sizeof(COLOR)*b->ninc); | 
| 689 | 
– | 
                                        /* compose direction transform */ | 
| 690 | 
– | 
        for (i = 3; i--; ) { | 
| 691 | 
– | 
                mymat[i][0] = u[i]; | 
| 692 | 
– | 
                mymat[i][1] = v[i]; | 
| 693 | 
– | 
                mymat[i][2] = w[i]; | 
| 694 | 
– | 
                mymat[i][3] = 0.; | 
| 695 | 
– | 
        } | 
| 696 | 
– | 
        mymat[3][0] = mymat[3][1] = mymat[3][2] = 0.; | 
| 697 | 
– | 
        mymat[3][3] = 1.; | 
| 698 | 
– | 
        if (xm != NULL) | 
| 699 | 
– | 
                multmat4(mymat, xm, mymat); | 
| 700 | 
– | 
        for (i = 3; i--; ) {            /* make sure it's normalized */ | 
| 701 | 
– | 
                wt = 1./sqrt(   mymat[0][i]*mymat[0][i] + | 
| 702 | 
– | 
                                mymat[1][i]*mymat[1][i] + | 
| 703 | 
– | 
                                mymat[2][i]*mymat[2][i] ); | 
| 704 | 
– | 
                for (j = 3; j--; ) | 
| 705 | 
– | 
                        mymat[j][i] *= wt; | 
| 706 | 
– | 
        } | 
| 707 | 
– | 
        if (!invmat4(inmat, mymat))     /* need inverse as well */ | 
| 708 | 
– | 
                error(INTERNAL, "cannot invert BSDF transform"); | 
| 709 | 
– | 
        newdist(nalt*nazi);             /* resample distribution */ | 
| 710 | 
– | 
        for (i = b->ninc; i--; ) { | 
| 711 | 
– | 
                int     direct_out = -1; | 
| 712 | 
– | 
                COLOR   cdir; | 
| 713 | 
– | 
                getBSDF_incvec(dv, b, i);       /* compute incident irrad. */ | 
| 714 | 
– | 
                multv3(dv, dv, mymat); | 
| 715 | 
– | 
                if (dv[2] < 0.0) { | 
| 716 | 
– | 
                        dv[0] = -dv[0]; dv[1] = -dv[1]; dv[2] = -dv[2]; | 
| 717 | 
– | 
                        direct_out += (direct_discount != NULL); | 
| 718 | 
– | 
                } | 
| 719 | 
– | 
                wt = getBSDF_incohm(b, i); | 
| 720 | 
– | 
                wt *= dv[2];                    /* solid_angle*cosine(theta) */ | 
| 721 | 
– | 
                cp = &idist[3*i]; | 
| 722 | 
– | 
                copycolor(cinc, cp); | 
| 723 | 
– | 
                scalecolor(cinc, wt); | 
| 724 | 
– | 
                if (!direct_out) {              /* discount direct contr. */ | 
| 725 | 
– | 
                        cp = &direct_discount[3*i]; | 
| 726 | 
– | 
                        copycolor(cdir, cp); | 
| 727 | 
– | 
                        scalecolor(cdir, -wt); | 
| 728 | 
– | 
                        direct_out = flatindex(dv, nalt, nazi); | 
| 729 | 
– | 
                } | 
| 730 | 
– | 
                for (k = nalt; k--; )           /* loop over distribution */ | 
| 731 | 
– | 
                  for (j = nazi; j--; ) { | 
| 732 | 
– | 
                    int rstart = random(); | 
| 733 | 
– | 
                    for (c = NBSDFSAMPS; c--; ) { | 
| 734 | 
– | 
                        double  sp[2]; | 
| 735 | 
– | 
                        multisamp(sp, 2, urand(rstart+c)); | 
| 736 | 
– | 
                        flatdir(dv, (k + sp[0])/nalt, | 
| 737 | 
– | 
                                        (j + .5 - sp[1])/nazi); | 
| 738 | 
– | 
                        multv3(dv, dv, inmat); | 
| 739 | 
– | 
                                                /* evaluate BSDF @ outgoing */ | 
| 740 | 
– | 
                        o = getBSDF_outndx(b, dv); | 
| 741 | 
– | 
                        if (o < 0) { | 
| 742 | 
– | 
                                nout++; | 
| 743 | 
– | 
                                continue; | 
| 744 | 
– | 
                        } | 
| 745 | 
– | 
                        wt = BSDF_value(b, i, o) * (1./NBSDFSAMPS); | 
| 746 | 
– | 
                        copycolor(col, cinc); | 
| 747 | 
– | 
                        o = k*nazi + j; | 
| 748 | 
– | 
                        if (o == direct_out) | 
| 749 | 
– | 
                                addcolor(col, cdir);    /* minus direct */ | 
| 750 | 
– | 
                        scalecolor(col, wt); | 
| 751 | 
– | 
                        cp = &distarr[3*o]; | 
| 752 | 
– | 
                        addcolor(cp, col);      /* sum into distribution */ | 
| 753 | 
– | 
                    } | 
| 754 | 
– | 
                  } | 
| 755 | 
– | 
        } | 
| 756 | 
– | 
        free(idist);                    /* free temp space */ | 
| 757 | 
– | 
        if (nout) { | 
| 758 | 
– | 
                sprintf(errmsg, "missing %.1f%% of BSDF directions", | 
| 759 | 
– | 
                                100.*nout/(b->ninc*nalt*nazi*NBSDFSAMPS)); | 
| 760 | 
– | 
                error(WARNING, errmsg); | 
| 761 | 
– | 
        } | 
| 538 | 
  | 
} |