| 252 | 
  | 
        return(NULL);   /* gratis return */ | 
| 253 | 
  | 
} | 
| 254 | 
  | 
 | 
| 255 | 
< | 
/* Multiply two matrices (or a matrix and a vector) and allocate the result*/ | 
| 255 | 
> | 
/* Scale a matrix by a single value */ | 
| 256 | 
  | 
static CMATRIX * | 
| 257 | 
+ | 
cm_scale(const CMATRIX *cm1, const COLOR sca) | 
| 258 | 
+ | 
{ | 
| 259 | 
+ | 
        CMATRIX *cmr; | 
| 260 | 
+ | 
        int     dr, dc; | 
| 261 | 
+ | 
 | 
| 262 | 
+ | 
        cmr = cm_alloc(cm1->nrows, cm1->ncols); | 
| 263 | 
+ | 
        if (cmr == NULL) | 
| 264 | 
+ | 
                return(NULL); | 
| 265 | 
+ | 
        for (dr = 0; dr < cmr->nrows; dr++) | 
| 266 | 
+ | 
            for (dc = 0; dc < cmr->ncols; dc++) { | 
| 267 | 
+ | 
                const COLORV    *sp = cm_lval(cm1,dr,dc); | 
| 268 | 
+ | 
                COLORV          *dp = cm_lval(cmr,dr,dc); | 
| 269 | 
+ | 
                dp[0] = sp[0] * sca[0]; | 
| 270 | 
+ | 
                dp[1] = sp[1] * sca[1]; | 
| 271 | 
+ | 
                dp[2] = sp[2] * sca[2]; | 
| 272 | 
+ | 
            } | 
| 273 | 
+ | 
        return(cmr); | 
| 274 | 
+ | 
} | 
| 275 | 
+ | 
 | 
| 276 | 
+ | 
/* Multiply two matrices (or a matrix and a vector) and allocate the result */ | 
| 277 | 
+ | 
static CMATRIX * | 
| 278 | 
  | 
cm_multiply(const CMATRIX *cm1, const CMATRIX *cm2) | 
| 279 | 
  | 
{ | 
| 280 | 
  | 
        CMATRIX *cmr; | 
| 322 | 
  | 
        int     nbadohm = 0; | 
| 323 | 
  | 
        int     nneg = 0; | 
| 324 | 
  | 
        int     r, c; | 
| 325 | 
< | 
                                        /* reciprocity is "transparent" */ | 
| 325 | 
> | 
                                        /* loop over incident angles */ | 
| 326 | 
  | 
        for (c = 0; c < cm->ncols; c++) { | 
| 327 | 
  | 
                const double    dom = mBSDF_incohm(bsdf,c); | 
| 328 | 
  | 
                                        /* projected solid angle */ | 
| 351 | 
  | 
        return(cm); | 
| 352 | 
  | 
} | 
| 353 | 
  | 
 | 
| 354 | 
+ | 
/* Convert between input and output indices for reciprocity */ | 
| 355 | 
+ | 
static int | 
| 356 | 
+ | 
recip_out_from_in(const SDMat *bsdf, int in_recip) | 
| 357 | 
+ | 
{ | 
| 358 | 
+ | 
        FVECT   v; | 
| 359 | 
+ | 
 | 
| 360 | 
+ | 
        if (!mBSDF_incvec(v, bsdf, in_recip+.5)) | 
| 361 | 
+ | 
                return(in_recip);               /* XXX should be error! */ | 
| 362 | 
+ | 
        v[2] = -v[2]; | 
| 363 | 
+ | 
        return(mBSDF_outndx(bsdf, v)); | 
| 364 | 
+ | 
} | 
| 365 | 
+ | 
 | 
| 366 | 
+ | 
/* Convert between output and input indices for reciprocity */ | 
| 367 | 
+ | 
static int | 
| 368 | 
+ | 
recip_in_from_out(const SDMat *bsdf, int out_recip) | 
| 369 | 
+ | 
{ | 
| 370 | 
+ | 
        FVECT   v; | 
| 371 | 
+ | 
 | 
| 372 | 
+ | 
        if (!mBSDF_outvec(v, bsdf, out_recip+.5)) | 
| 373 | 
+ | 
                return(out_recip);              /* XXX should be error! */ | 
| 374 | 
+ | 
        v[2] = -v[2]; | 
| 375 | 
+ | 
        return(mBSDF_incndx(bsdf, v)); | 
| 376 | 
+ | 
} | 
| 377 | 
+ | 
 | 
| 378 | 
+ | 
/* Convert a BSDF to our matrix representation, applying reciprocity */ | 
| 379 | 
+ | 
static CMATRIX * | 
| 380 | 
+ | 
cm_bsdf_recip(const COLOR bsdfLamb, const COLOR specCol, const SDMat *bsdf) | 
| 381 | 
+ | 
{ | 
| 382 | 
+ | 
        CMATRIX *cm = cm_alloc(bsdf->ninc, bsdf->nout); | 
| 383 | 
+ | 
        int     nbadohm = 0; | 
| 384 | 
+ | 
        int     nneg = 0; | 
| 385 | 
+ | 
        int     r, c; | 
| 386 | 
+ | 
                                        /* loop over incident angles */ | 
| 387 | 
+ | 
        for (c = 0; c < cm->ncols; c++) { | 
| 388 | 
+ | 
                const int       ro = recip_out_from_in(bsdf,c); | 
| 389 | 
+ | 
                const double    dom = mBSDF_outohm(bsdf,ro); | 
| 390 | 
+ | 
                                        /* projected solid angle */ | 
| 391 | 
+ | 
                nbadohm += (dom <= 0); | 
| 392 | 
+ | 
 | 
| 393 | 
+ | 
                for (r = 0; r < cm->nrows; r++) { | 
| 394 | 
+ | 
                        const int       ri = recip_in_from_out(bsdf,r); | 
| 395 | 
+ | 
                        float           f = mBSDF_value(bsdf,ri,ro); | 
| 396 | 
+ | 
                        COLORV          *mp = cm_lval(cm,r,c); | 
| 397 | 
+ | 
                                        /* check BSDF value */ | 
| 398 | 
+ | 
                        if ((f <= 0) | (dom <= 0)) { | 
| 399 | 
+ | 
                                nneg += (f < -FTINY); | 
| 400 | 
+ | 
                                f = .0f; | 
| 401 | 
+ | 
                        } | 
| 402 | 
+ | 
                        copycolor(mp, specCol); | 
| 403 | 
+ | 
                        scalecolor(mp, f); | 
| 404 | 
+ | 
                        addcolor(mp, bsdfLamb); | 
| 405 | 
+ | 
                        scalecolor(mp, dom); | 
| 406 | 
+ | 
                } | 
| 407 | 
+ | 
        } | 
| 408 | 
+ | 
        if (nneg | nbadohm) { | 
| 409 | 
+ | 
                sprintf(errmsg, | 
| 410 | 
+ | 
                    "BTDF has %d negatives and %d bad incoming solid angles", | 
| 411 | 
+ | 
                                nneg, nbadohm); | 
| 412 | 
+ | 
                error(WARNING, errmsg); | 
| 413 | 
+ | 
        } | 
| 414 | 
+ | 
        return(cm); | 
| 415 | 
+ | 
} | 
| 416 | 
+ | 
 | 
| 417 | 
  | 
/* Load and convert a matrix BSDF from the given XML file */ | 
| 418 | 
  | 
static CMATRIX * | 
| 419 | 
< | 
cm_loadBSDF(char *fname) | 
| 419 | 
> | 
cm_loadBSDF(char *fname, COLOR cLamb) | 
| 420 | 
  | 
{ | 
| 421 | 
< | 
        CMATRIX *Tmat; | 
| 422 | 
< | 
        char    *fpath; | 
| 423 | 
< | 
        SDError ec; | 
| 424 | 
< | 
        SDData  myBSDF; | 
| 425 | 
< | 
        COLOR   bsdfLamb, specCol; | 
| 421 | 
> | 
        CMATRIX         *Tmat; | 
| 422 | 
> | 
        char            *fpath; | 
| 423 | 
> | 
        int             recip; | 
| 424 | 
> | 
        SDError         ec; | 
| 425 | 
> | 
        SDData          myBSDF; | 
| 426 | 
> | 
        SDSpectralDF    *tdf; | 
| 427 | 
> | 
        COLOR           bsdfLamb, specCol; | 
| 428 | 
  | 
                                        /* find path to BSDF file */ | 
| 429 | 
  | 
        fpath = getpath(fname, getrlibpath(), R_OK); | 
| 430 | 
  | 
        if (fpath == NULL) { | 
| 435 | 
  | 
        ec = SDloadFile(&myBSDF, fpath); | 
| 436 | 
  | 
        if (ec) | 
| 437 | 
  | 
                error(USER, transSDError(ec)); | 
| 438 | 
< | 
        if (myBSDF.tf == NULL || myBSDF.tf->ncomp != 1 || | 
| 439 | 
< | 
                        myBSDF.tf->comp[0].func != &SDhandleMtx) { | 
| 438 | 
> | 
        ccy2rgb(&myBSDF.tLamb.spec, myBSDF.tLamb.cieY/PI, bsdfLamb); | 
| 439 | 
> | 
        recip = (myBSDF.tf == NULL); | 
| 440 | 
> | 
        tdf = recip ? myBSDF.tb : myBSDF.tf; | 
| 441 | 
> | 
        if (tdf == NULL) {              /* no non-Lambertian transmission? */ | 
| 442 | 
> | 
                if (cLamb != NULL) | 
| 443 | 
> | 
                        copycolor(cLamb, bsdfLamb); | 
| 444 | 
> | 
                SDfreeBSDF(&myBSDF); | 
| 445 | 
> | 
                return(NULL); | 
| 446 | 
> | 
        } | 
| 447 | 
> | 
        if (tdf->ncomp != 1 || tdf->comp[0].func != &SDhandleMtx) { | 
| 448 | 
  | 
                sprintf(errmsg, "unsupported BSDF '%s'", fpath); | 
| 449 | 
  | 
                error(USER, errmsg); | 
| 450 | 
  | 
        } | 
| 451 | 
  | 
                                        /* convert BTDF to matrix */ | 
| 452 | 
< | 
        ccy2rgb(&myBSDF.tLamb.spec, myBSDF.tLamb.cieY/PI, bsdfLamb); | 
| 453 | 
< | 
        ccy2rgb(&myBSDF.tf->comp[0].cspec[0], 1., specCol); | 
| 454 | 
< | 
        Tmat = cm_bsdf(bsdfLamb, specCol, (SDMat *)myBSDF.tf->comp[0].dist); | 
| 452 | 
> | 
        ccy2rgb(&tdf->comp[0].cspec[0], 1., specCol); | 
| 453 | 
> | 
        Tmat = recip ? cm_bsdf_recip(bsdfLamb, specCol, (SDMat *)tdf->comp[0].dist) | 
| 454 | 
> | 
                        : cm_bsdf(bsdfLamb, specCol, (SDMat *)tdf->comp[0].dist); | 
| 455 | 
> | 
        if (cLamb != NULL)              /* Lambertian is included */ | 
| 456 | 
> | 
                setcolor(cLamb, .0, .0, .0); | 
| 457 | 
  | 
                                        /* free BSDF and return */ | 
| 458 | 
  | 
        SDfreeBSDF(&myBSDF); | 
| 459 | 
  | 
        return(Tmat); | 
| 478 | 
  | 
                int             dt, xr, yr; | 
| 479 | 
  | 
                COLORV          *psp; | 
| 480 | 
  | 
                                                        /* check for zero */ | 
| 481 | 
< | 
                if ((scv[RED] == 0) & (scv[GRN] == 0) & (scv[BLU] == 0)) | 
| 481 | 
> | 
                if ((scv[RED] == 0) & (scv[GRN] == 0) & (scv[BLU] == 0) && | 
| 482 | 
> | 
                                (myDT != DTfromHeader) | (i < cv->nrows-1)) | 
| 483 | 
  | 
                        continue; | 
| 484 | 
  | 
                                                        /* open next picture */ | 
| 485 | 
  | 
                sprintf(fname, fspec, i); | 
| 562 | 
  | 
int | 
| 563 | 
  | 
main(int argc, char *argv[]) | 
| 564 | 
  | 
{ | 
| 565 | 
< | 
        CMATRIX                 *cvec; | 
| 565 | 
> | 
        CMATRIX                 *cvec;          /* component vector result */ | 
| 566 | 
  | 
 | 
| 567 | 
  | 
        progname = argv[0]; | 
| 568 | 
  | 
 | 
| 575 | 
  | 
 | 
| 576 | 
  | 
        if (argc > 3) {                         /* VTDs expression */ | 
| 577 | 
  | 
                CMATRIX *svec, *Dmat, *Tmat, *ivec; | 
| 578 | 
+ | 
                COLOR   tLamb; | 
| 579 | 
  | 
                                                /* get sky vector */ | 
| 580 | 
  | 
                svec = cm_load(argv[4], 0, 1, DTascii); | 
| 581 | 
  | 
                                                /* load BSDF */ | 
| 582 | 
< | 
                Tmat = cm_loadBSDF(argv[2]); | 
| 582 | 
> | 
                Tmat = cm_loadBSDF(argv[2], tLamb); | 
| 583 | 
  | 
                                                /* load Daylight matrix */ | 
| 584 | 
< | 
                Dmat = cm_load(argv[3], Tmat->ncols, svec->nrows, DTfromHeader); | 
| 584 | 
> | 
                Dmat = cm_load(argv[3], Tmat==NULL ? 0 : Tmat->ncols, | 
| 585 | 
> | 
                                        svec->nrows, DTfromHeader); | 
| 586 | 
  | 
                                                /* multiply vector through */ | 
| 587 | 
  | 
                ivec = cm_multiply(Dmat, svec); | 
| 588 | 
  | 
                cm_free(Dmat); cm_free(svec); | 
| 589 | 
< | 
                cvec = cm_multiply(Tmat, ivec); /* cvec = component vector */ | 
| 590 | 
< | 
                cm_free(Tmat); cm_free(ivec); | 
| 589 | 
> | 
                if (Tmat == NULL) {             /* diffuse only */ | 
| 590 | 
> | 
                        cvec = cm_scale(ivec, tLamb); | 
| 591 | 
> | 
                } else {                        /* else apply BTDF matrix */ | 
| 592 | 
> | 
                        cvec = cm_multiply(Tmat, ivec); | 
| 593 | 
> | 
                        cm_free(Tmat);  | 
| 594 | 
> | 
                } | 
| 595 | 
> | 
                cm_free(ivec); | 
| 596 | 
  | 
        } else {                                /* else just use sky vector */ | 
| 597 | 
  | 
                cvec = cm_load(argv[2], 0, 1, DTascii); | 
| 598 | 
  | 
        } |