ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/cv/bsdf2rad.c
(Generate patch)

Comparing ray/src/cv/bsdf2rad.c (file contents):
Revision 2.13 by greg, Mon Sep 15 00:54:03 2014 UTC vs.
Revision 2.35 by greg, Thu Aug 2 18:33:42 2018 UTC

# Line 5 | Line 5 | static const char RCSid[] = "$Id$";
5   *  Plot 3-D BSDF output based on scattering interpolant or XML representation
6   */
7  
8 #define _USE_MATH_DEFINES
8   #include <stdio.h>
9   #include <string.h>
10   #include <stdlib.h>
11 < #include <math.h>
12 < #include "rtprocess.h"
11 > #include "paths.h"
12 > #include "rtmath.h"
13 > #include "resolu.h"
14   #include "bsdfrep.h"
15  
16 < const float     colarr[6][3] = {
17 <                .7, 1., .7,
18 <                1., .7, .7,
19 <                .7, .7, 1.,
20 <                1., .5, 1.,
21 <                1., 1., .5,
22 <                .5, 1., 1.
23 <        };
16 > #ifndef NINCIDENT
17 > #define NINCIDENT       37              /* number of samples/hemisphere */
18 > #endif
19 > #ifndef GRIDSTEP
20 > #define GRIDSTEP        2               /* our grid step size */
21 > #endif
22 > #define SAMPRES         (GRIDRES/GRIDSTEP)
23  
24 < #ifdef _WIN32
25 < char    validf[] = "-e \"valid(s,t)=X`SYS(s,t)^2+Y`SYS(s,t)^2+Z`SYS(s,t)^2-1e-7\"";
24 > int     front_comp = 0;                 /* front component flags (SDsamp*) */
25 > int     back_comp = 0;                  /* back component flags */
26 > double  overall_min = 1./PI;            /* overall minimum BSDF value */
27 > double  min_log10;                      /* smallest log10 value for plotting */
28 > double  overall_max = .0;               /* overall maximum BSDF value */
29 >
30 > char    ourTempDir[TEMPLEN+1] = "";     /* our temporary directory */
31 >
32 > const char      frpref[] = "rf";
33 > const char      ftpref[] = "tf";
34 > const char      brpref[] = "rb";
35 > const char      btpref[] = "tb";
36 > const char      dsuffix[] = ".txt";
37 >
38 > const char      sph_fmat[] = "fBSDFmat";
39 > const char      sph_bmat[] = "bBSDFmat";
40 > const double    sph_rad = 10.;
41 > const double    sph_xoffset = 15.;
42 >
43 > #define bsdf_rad        (sph_rad*.25)
44 > #define arrow_rad       (bsdf_rad*.015)
45 >
46 > #define FEQ(a,b)        ((a)-(b) <= 1e-7 && (b)-(a) <= 1e-7)
47 >
48 > #define set_minlog()    overall_min = (overall_min < 1e-5) ? 1e-5 : overall_min; \
49 >                                min_log10 = log10(overall_min) - .1
50 >
51 > char    *progname;
52 >
53 > /* Get Fibonacci sphere vector (0 to NINCIDENT-1) */
54 > static RREAL *
55 > get_ivector(FVECT iv, int i)
56 > {
57 >        const double    phistep = PI*(3. - 2.236067978);
58 >        double          r;
59 >
60 >        iv[2] = 1. - (i+.5)*(1./NINCIDENT);
61 >        r = sqrt(1. - iv[2]*iv[2]);
62 >        iv[0] = r * cos((i+1.)*phistep);
63 >        iv[1] = r * sin((i+1.)*phistep);
64 >
65 >        return(iv);
66 > }
67 >
68 > /* Convert incident vector into sphere position */
69 > static RREAL *
70 > cvt_sposition(FVECT sp, const FVECT iv, int inc_side)
71 > {
72 >        sp[0] = -iv[0]*sph_rad + inc_side*sph_xoffset;
73 >        sp[1] = -iv[1]*sph_rad;
74 >        sp[2] = iv[2]*sph_rad;
75 >
76 >        return(sp);
77 > }
78 >
79 > /* Get temporary file name */
80 > static char *
81 > tfile_name(const char *prefix, const char *suffix, int i)
82 > {
83 >        static char     buf[256];
84 >
85 >        if (!ourTempDir[0]) {           /* create temporary directory */
86 >                mktemp(strcpy(ourTempDir,TEMPLATE));
87 >                if (mkdir(ourTempDir, 0777) < 0) {
88 >                        perror("mkdir");
89 >                        exit(1);
90 >                }
91 >        }
92 >        if (!prefix) prefix = "T";
93 >        if (!suffix) suffix = "";
94 >        sprintf(buf, "%s/%s%03d%s", ourTempDir, prefix, i, suffix);
95 >        return(buf);
96 > }
97 >
98 > /* Remove temporary directory & contents */
99 > static void
100 > cleanup_tmp(void)
101 > {
102 >        char    buf[128];
103 >
104 >        if (!ourTempDir[0])
105 >                return;
106 > #if defined(_WIN32) || defined(_WIN64)
107 >        sprintf(buf, "RMDIR %s /S /Q", ourTempDir);
108   #else
109 < char    validf[] = "-e 'valid(s,t)=X`SYS(s,t)^2+Y`SYS(s,t)^2+Z`SYS(s,t)^2-1e-7'";
109 >        sprintf(buf, "rm -rf %s", ourTempDir);
110   #endif
111 +        system(buf);
112 + }
113  
114 < char    *progname;
114 > /* Run the specified command, returning 1 if OK */
115 > static int
116 > run_cmd(const char *cmd)
117 > {
118 >        fflush(stdout);
119 >        if (system(cmd)) {
120 >                fprintf(stderr, "%s: error running: %s\n", progname, cmd);
121 >                return(0);
122 >        }
123 >        return(1);
124 > }
125  
126 < /* Produce a Radiance model plotting the indicated incident direction(s) */
126 > /* Plot surface points for the given BSDF incident angle */
127 > static int
128 > plotBSDF(const char *fname, const FVECT ivec, int dfl, const SDData *sd)
129 > {
130 >        FILE            *fp = fopen(fname, "w");
131 >        int             i, j;
132 >
133 >        if (fp == NULL) {
134 >                fprintf(stderr, "%s: cannot open '%s' for writing\n",
135 >                                progname, fname);
136 >                return(0);
137 >        }
138 >        if (ivec[2] > 0) {
139 >                input_orient = 1;
140 >                output_orient = dfl&SDsampR ? 1 : -1;
141 >        } else {
142 >                input_orient = -1;
143 >                output_orient = dfl&SDsampR ? -1 : 1;
144 >        }
145 >        for (i = SAMPRES; i--; )
146 >            for (j = 0; j < SAMPRES; j++) {
147 >                FVECT   ovec;
148 >                SDValue sval;
149 >                double  bsdf;
150 >                ovec_from_pos(ovec, i*GRIDSTEP, j*GRIDSTEP);
151 >                if (SDreportError(SDevalBSDF(&sval, ovec,
152 >                                                ivec, sd), stderr))
153 >                        return(0);
154 >                if (sval.cieY > overall_max)
155 >                        overall_max = sval.cieY;
156 >                bsdf = (sval.cieY < overall_min) ? overall_min : sval.cieY;
157 >                bsdf = log10(bsdf) - min_log10;
158 >                fprintf(fp, "%.5f %.5f %.5f\n",
159 >                                ovec[0]*bsdf, ovec[1]*bsdf, ovec[2]*bsdf);
160 >            }
161 >        if (fclose(fp) == EOF) {
162 >                fprintf(stderr, "%s: error writing data to '%s'\n",
163 >                                progname, fname);
164 >                return(0);
165 >        }
166 >        return(1);
167 > }
168 >
169 > /* Build BSDF values from loaded XML file */
170 > static int
171 > build_wBSDF(const SDData *sd)
172 > {
173 >        FVECT   ivec;
174 >        int     i;
175 >
176 >        if (front_comp & SDsampR)
177 >                for (i = 0; i < NINCIDENT; i++) {
178 >                        get_ivector(ivec, i);
179 >                        if (!plotBSDF(tfile_name(frpref, dsuffix, i),
180 >                                        ivec, SDsampR, sd))
181 >                                return(0);
182 >                }
183 >        if (front_comp & SDsampT)
184 >                for (i = 0; i < NINCIDENT; i++) {
185 >                        get_ivector(ivec, i);
186 >                        if (!plotBSDF(tfile_name(ftpref, dsuffix, i),
187 >                                        ivec, SDsampT, sd))
188 >                                return(0);
189 >                }
190 >        if (back_comp & SDsampR)
191 >                for (i = 0; i < NINCIDENT; i++) {
192 >                        get_ivector(ivec, i);
193 >                        ivec[0] = -ivec[0]; ivec[2] = -ivec[2];
194 >                        if (!plotBSDF(tfile_name(brpref, dsuffix, i),
195 >                                        ivec, SDsampR, sd))
196 >                                return(0);
197 >                }
198 >        if (back_comp & SDsampT)
199 >                for (i = 0; i < NINCIDENT; i++) {
200 >                        get_ivector(ivec, i);
201 >                        ivec[0] = -ivec[0]; ivec[2] = -ivec[2];
202 >                        if (!plotBSDF(tfile_name(btpref, dsuffix, i),
203 >                                        ivec, SDsampT, sd))
204 >                                return(0);
205 >                }
206 >        return(1);
207 > }
208 >
209 > /* Plot surface points using radial basis function */
210 > static int
211 > plotRBF(const char *fname, const RBFNODE *rbf)
212 > {
213 >        FILE            *fp = fopen(fname, "w");
214 >        int             i, j;
215 >
216 >        if (fp == NULL) {
217 >                fprintf(stderr, "%s: cannot open '%s' for writing\n",
218 >                                progname, fname);
219 >                return(0);
220 >        }
221 >        for (i = SAMPRES; i--; )
222 >            for (j = 0; j < SAMPRES; j++) {
223 >                FVECT   ovec;
224 >                double  bsdf;
225 >                ovec_from_pos(ovec, i*GRIDSTEP, j*GRIDSTEP);
226 >                bsdf = eval_rbfrep(rbf, ovec);
227 >                if (bsdf > overall_max)
228 >                        overall_max = bsdf;
229 >                else if (bsdf < overall_min)
230 >                        bsdf = overall_min;
231 >                bsdf = log10(bsdf) - min_log10;
232 >                fprintf(fp, "%.5f %.5f %.5f\n",
233 >                                ovec[0]*bsdf, ovec[1]*bsdf, ovec[2]*bsdf);
234 >            }
235 >        if (fclose(fp) == EOF) {
236 >                fprintf(stderr, "%s: error writing data to '%s'\n",
237 >                                progname, fname);
238 >                return(0);
239 >        }
240 >        return(1);
241 > }
242 >
243 > /* Build BSDF values from scattering interpolant representation */
244 > static int
245 > build_wRBF(void)
246 > {
247 >        const char      *pref;
248 >        int             i;
249 >
250 >        if (input_orient > 0) {
251 >                if (output_orient > 0)
252 >                        pref = frpref;
253 >                else
254 >                        pref = ftpref;
255 >        } else if (output_orient < 0)
256 >                pref = brpref;
257 >        else
258 >                pref = btpref;
259 >
260 >        for (i = 0; i < NINCIDENT; i++) {
261 >                FVECT   ivec;
262 >                RBFNODE *rbf;
263 >                get_ivector(ivec, i);
264 >                if (input_orient < 0) {
265 >                        ivec[0] = -ivec[0]; ivec[2] = -ivec[2];
266 >                }
267 >                rbf = advect_rbf(ivec, 15000);
268 >                if (!plotRBF(tfile_name(pref, dsuffix, i), rbf))
269 >                        return(0);
270 >                if (rbf) free(rbf);
271 >        }
272 >        return(1);                              /* next call frees */
273 > }
274 >
275 > /* Put out mirror arrow for the given incident vector */
276 > static void
277 > put_mirror_arrow(const FVECT origin, const FVECT nrm)
278 > {
279 >        const double    arrow_len = 1.2*bsdf_rad;
280 >        const double    tip_len = 0.2*bsdf_rad;
281 >        static int      cnt = 1;
282 >        FVECT           refl;
283 >        int             i;
284 >
285 >        refl[0] = 2.*nrm[2]*nrm[0];
286 >        refl[1] = 2.*nrm[2]*nrm[1];
287 >        refl[2] = 2.*nrm[2]*nrm[2] - 1.;
288 >
289 >        printf("\n# Mirror arrow #%d\n", cnt);
290 >        printf("\nshaft_mat cylinder inc_dir%d\n0\n0\n7", cnt);
291 >        printf("\n\t%f %f %f\n\t%f %f %f\n\t%f\n",
292 >                        origin[0], origin[1], origin[2]+arrow_len,
293 >                        origin[0], origin[1], origin[2],
294 >                        arrow_rad);
295 >        printf("\nshaft_mat cylinder mir_dir%d\n0\n0\n7", cnt);
296 >        printf("\n\t%f %f %f\n\t%f %f %f\n\t%f\n",
297 >                        origin[0], origin[1], origin[2],
298 >                        origin[0] + arrow_len*refl[0],
299 >                        origin[1] + arrow_len*refl[1],
300 >                        origin[2] + arrow_len*refl[2],
301 >                        arrow_rad);
302 >        printf("\ntip_mat cone mir_tip%d\n0\n0\n8", cnt);
303 >        printf("\n\t%f %f %f\n\t%f %f %f\n\t%f 0\n",
304 >                        origin[0] + (arrow_len-.5*tip_len)*refl[0],
305 >                        origin[1] + (arrow_len-.5*tip_len)*refl[1],
306 >                        origin[2] + (arrow_len-.5*tip_len)*refl[2],
307 >                        origin[0] + (arrow_len+.5*tip_len)*refl[0],
308 >                        origin[1] + (arrow_len+.5*tip_len)*refl[1],
309 >                        origin[2] + (arrow_len+.5*tip_len)*refl[2],
310 >                        2.*arrow_rad);
311 >        ++cnt;
312 > }
313 >
314 > /* Put out transmitted direction arrow for the given incident vector */
315 > static void
316 > put_trans_arrow(const FVECT origin)
317 > {
318 >        const double    arrow_len = 1.2*bsdf_rad;
319 >        const double    tip_len = 0.2*bsdf_rad;
320 >        static int      cnt = 1;
321 >        int             i;
322 >
323 >        printf("\n# Transmission arrow #%d\n", cnt);
324 >        printf("\nshaft_mat cylinder trans_dir%d\n0\n0\n7", cnt);
325 >        printf("\n\t%f %f %f\n\t%f %f %f\n\t%f\n",
326 >                        origin[0], origin[1], origin[2],
327 >                        origin[0], origin[1], origin[2]-arrow_len,
328 >                        arrow_rad);
329 >        printf("\ntip_mat cone trans_tip%d\n0\n0\n8", cnt);
330 >        printf("\n\t%f %f %f\n\t%f %f %f\n\t%f 0\n",
331 >                        origin[0], origin[1], origin[2]-arrow_len+.5*tip_len,
332 >                        origin[0], origin[1], origin[2]-arrow_len-.5*tip_len,
333 >                        2.*arrow_rad);
334 >        ++cnt;
335 > }
336 >
337 > /* Compute rotation (x,y,z) => (xp,yp,zp) */
338 > static int
339 > addrot(char *xf, const FVECT xp, const FVECT yp, const FVECT zp)
340 > {
341 >        int     n = 0;
342 >        double  theta;
343 >
344 >        if (yp[2]*yp[2] + zp[2]*zp[2] < 2.*FTINY*FTINY) {
345 >                /* Special case for X' along Z-axis */
346 >                theta = -atan2(yp[0], yp[1]);
347 >                sprintf(xf, " -ry %f -rz %f",
348 >                                xp[2] < 0.0 ? 90.0 : -90.0,
349 >                                theta*(180./PI));
350 >                return(4);
351 >        }
352 >        theta = atan2(yp[2], zp[2]);
353 >        if (!FEQ(theta,0.0)) {
354 >                sprintf(xf, " -rx %f", theta*(180./PI));
355 >                while (*xf) ++xf;
356 >                n += 2;
357 >        }
358 >        theta = Asin(-xp[2]);
359 >        if (!FEQ(theta,0.0)) {
360 >                sprintf(xf, " -ry %f", theta*(180./PI));
361 >                while (*xf) ++xf;
362 >                n += 2;
363 >        }
364 >        theta = atan2(xp[1], xp[0]);
365 >        if (!FEQ(theta,0.0)) {
366 >                sprintf(xf, " -rz %f", theta*(180./PI));
367 >                /* while (*xf) ++xf; */
368 >                n += 2;
369 >        }
370 >        return(n);
371 > }
372 >
373 > /* Put out BSDF surfaces */
374 > static int
375 > put_BSDFs(void)
376 > {
377 >        const double    scalef = bsdf_rad/(log10(overall_max) - min_log10);
378 >        FVECT           ivec, sorg, nrm, upv;
379 >        RREAL           vMtx[3][3];
380 >        char            *fname;
381 >        char            cmdbuf[256];
382 >        char            rotargs[64];
383 >        int             nrota;
384 >        int             i;
385 >
386 >        printf("\n# Gensurf output corresponding to %d incident directions\n",
387 >                        NINCIDENT);
388 >
389 >        printf("\nvoid glow tip_mat\n0\n0\n4 1 0 1 0\n");
390 >        printf("\nvoid mixfunc shaft_mat\n4 tip_mat void 0.25 .\n0\n0\n");
391 >
392 >        for (i = 0; i < NINCIDENT; i++) {
393 >                get_ivector(ivec, i);
394 >                nrm[0] = -ivec[0]; nrm[1] = -ivec[1]; nrm[2] = ivec[2];
395 >                upv[0] = nrm[0]*nrm[1]*(nrm[2] - 1.);
396 >                upv[1] = nrm[0]*nrm[0] + nrm[1]*nrm[1]*nrm[2];
397 >                upv[2] = -nrm[1]*(nrm[0]*nrm[0] + nrm[1]*nrm[1]);
398 >                if (SDcompXform(vMtx, nrm, upv) != SDEnone)
399 >                        continue;
400 >                nrota = addrot(rotargs, vMtx[0], vMtx[1], vMtx[2]);
401 >                if (front_comp) {
402 >                        cvt_sposition(sorg, ivec, 1);
403 >                        printf("\nvoid colorfunc scale_pat\n");
404 >                        printf("10 bsdf_red bsdf_grn bsdf_blu bsdf2rad.cal\n");
405 >                        printf("\t-s %f -t %f %f %f\n0\n0\n",
406 >                                        bsdf_rad, sorg[0], sorg[1], sorg[2]);
407 >                        printf("\nscale_pat glow scale_mat\n0\n0\n4 1 1 1 0\n");
408 >                }
409 >                if (front_comp & SDsampR) {
410 >                        put_mirror_arrow(sorg, nrm);
411 >                        fname = tfile_name(frpref, dsuffix, i);
412 >                        sprintf(cmdbuf,
413 >                "gensurf scale_mat %s%d %s %s %s %d %d | xform %s -s %f -t %f %f %f",
414 >                                        frpref, i, fname, fname, fname, SAMPRES-1, SAMPRES-1,
415 >                                        rotargs, scalef, sorg[0], sorg[1], sorg[2]);
416 >                        if (!run_cmd(cmdbuf))
417 >                                return(0);
418 >                }
419 >                if (front_comp & SDsampT) {
420 >                        put_trans_arrow(sorg);
421 >                        fname = tfile_name(ftpref, dsuffix, i);
422 >                        sprintf(cmdbuf,
423 >                "gensurf scale_mat %s%d %s %s %s %d %d | xform -I %s -s %f -t %f %f %f",
424 >                                        ftpref, i, fname, fname, fname, SAMPRES-1, SAMPRES-1,
425 >                                        rotargs, scalef, sorg[0], sorg[1], sorg[2]);
426 >                        if (!run_cmd(cmdbuf))
427 >                                return(0);
428 >                }
429 >                if (back_comp) {
430 >                        cvt_sposition(sorg, ivec, -1);
431 >                        printf("\nvoid colorfunc scale_pat\n");
432 >                        printf("10 bsdf_red bsdf_grn bsdf_blu bsdf2rad.cal\n");
433 >                        printf("\t-s %f -t %f %f %f\n0\n0\n",
434 >                                        bsdf_rad, sorg[0], sorg[1], sorg[2]);
435 >                        printf("\nscale_pat glow scale_mat\n0\n0\n4 1 1 1 0\n");
436 >                }
437 >                if (back_comp & SDsampR) {
438 >                        put_mirror_arrow(sorg, nrm);
439 >                        fname = tfile_name(brpref, dsuffix, i);
440 >                        sprintf(cmdbuf,
441 >                "gensurf scale_mat %s%d %s %s %s %d %d | xform -I -ry 180 %s -s %f -t %f %f %f",
442 >                                        brpref, i, fname, fname, fname, SAMPRES-1, SAMPRES-1,
443 >                                        rotargs, scalef, sorg[0], sorg[1], sorg[2]);
444 >                        if (!run_cmd(cmdbuf))
445 >                                return(0);
446 >                }
447 >                if (back_comp & SDsampT) {
448 >                        put_trans_arrow(sorg);
449 >                        fname = tfile_name(btpref, dsuffix, i);
450 >                        sprintf(cmdbuf,
451 >                "gensurf scale_mat %s%d %s %s %s %d %d | xform -ry 180 %s -s %f -t %f %f %f",
452 >                                        btpref, i, fname, fname, fname, SAMPRES-1, SAMPRES-1,
453 >                                        rotargs, scalef, sorg[0], sorg[1], sorg[2]);
454 >                        if (!run_cmd(cmdbuf))
455 >                                return(0);
456 >                }
457 >        }
458 >        return(1);
459 > }
460 >
461 > /* Put our hemisphere material */
462 > static void
463 > put_matBSDF(const char *XMLfile)
464 > {
465 >        const char      *curdir = "./";
466 >
467 >        if (!XMLfile) {                 /* simple material */
468 >                printf("\n# Simplified material because we have no XML input\n");
469 >                printf("\nvoid brightfunc latlong\n2 latlong bsdf2rad.cal\n0\n0\n");
470 >                if ((front_comp|back_comp) & SDsampT)
471 >                        printf("\nlatlong trans %s\n0\n0\n7 .75 .75 .75 0 .04 .5 .8\n",
472 >                                        sph_fmat);
473 >                else
474 >                        printf("\nlatlong plastic %s\n0\n0\n5 .5 .5 .5 0 0\n",
475 >                                        sph_fmat);
476 >                printf("\ninherit alias %s %s\n", sph_bmat, sph_fmat);
477 >                return;
478 >        }
479 >        switch (XMLfile[0]) {           /* avoid RAYPATH search */
480 >        case '.':
481 >        case '~':
482 >        CASEDIRSEP:
483 >                curdir = "";
484 >                break;
485 >        case '\0':
486 >                fprintf(stderr, "%s: empty file name in put_matBSDF\n", progname);
487 >                exit(1);
488 >                break;
489 >        }
490 >        printf("\n# Actual BSDF materials for rendering the hemispheres\n");
491 >        printf("\nvoid BSDF BSDF_f\n6 0 \"%s%s\" upx upy upz bsdf2rad.cal\n0\n0\n",
492 >                        curdir, XMLfile);
493 >        printf("\nvoid plastic black\n0\n0\n5 0 0 0 0 0\n");
494 >        printf("\nvoid mixfunc %s\n4 BSDF_f black latlong bsdf2rad.cal\n0\n0\n",
495 >                        sph_fmat);
496 >        printf("\nvoid BSDF BSDF_b\n8 0 \"%s%s\" upx upy upz bsdf2rad.cal -ry 180\n0\n0\n",
497 >                        curdir, XMLfile);
498 >        printf("\nvoid mixfunc %s\n4 BSDF_b black latlong bsdf2rad.cal\n0\n0\n",
499 >                        sph_bmat);
500 > }
501 >
502 > /* Put out overhead parallel light source */
503 > static void
504 > put_source(void)
505 > {
506 >        printf("\n# Overhead parallel light source\n");
507 >        printf("\nvoid light bright\n0\n0\n3 2500 2500 2500\n");
508 >        printf("\nbright source light\n0\n0\n4 0 0 1 2\n");
509 >        printf("\n# Material used for labels\n");
510 >        printf("\nvoid trans vellum\n0\n0\n7 1 1 1 0 0 .5 0\n");
511 > }
512 >
513 > /* Put out hemisphere(s) */
514 > static void
515 > put_hemispheres(void)
516 > {
517 >        const int       nsegs = 131;
518 >
519 >        printf("\n# Hemisphere(s) for showing BSDF appearance (if XML file)\n");
520 >        if (front_comp) {
521 >                printf(
522 > "\n!genrev %s Front \"R*sin(A*t)\" \"R*cos(A*t)\" %d -e \"R:%g;A:%f\" -s | xform -t %g 0 0\n",
523 >                                sph_fmat, nsegs, sph_rad, 0.5*PI, sph_xoffset);
524 >                printf("\nvoid brighttext front_text\n3 helvet.fnt . FRONT\n0\n");
525 >                printf("12\n\t%f %f 0\n\t%f 0 0\n\t0 %f 0\n\t.01 1 -.1\n",
526 >                                -.22*sph_rad + sph_xoffset, -1.4*sph_rad,
527 >                                .35/5.*sph_rad, -1.6*.35/5.*sph_rad);
528 >                printf("\nfront_text alias front_label_mat vellum\n");
529 >                printf("\nfront_label_mat polygon front_label\n0\n0\n12");
530 >                printf("\n\t%f %f 0\n\t%f %f 0\n\t%f %f 0\n\t%f %f 0\n",
531 >                                -.25*sph_rad + sph_xoffset, -1.3*sph_rad,
532 >                                -.25*sph_rad + sph_xoffset, (-1.4-1.6*.35/5.-.1)*sph_rad,
533 >                                .25*sph_rad + sph_xoffset, (-1.4-1.6*.35/5.-.1)*sph_rad,
534 >                                .25*sph_rad + sph_xoffset, -1.3*sph_rad );
535 >        }
536 >        if (back_comp) {
537 >                printf(
538 > "\n!genrev %s Back \"R*cos(A*t)\" \"R*sin(A*t)\" %d -e \"R:%g;A:%f\" -s | xform -t %g 0 0\n",
539 >                                sph_bmat, nsegs, sph_rad, 0.5*PI, -sph_xoffset);
540 >                printf("\nvoid brighttext back_text\n3 helvet.fnt . BACK\n0\n");
541 >                printf("12\n\t%f %f 0\n\t%f 0 0\n\t0 %f 0\n\t.01 1 -.1\n",
542 >                                -.22*sph_rad - sph_xoffset, -1.4*sph_rad,
543 >                                .35/4.*sph_rad, -1.6*.35/4.*sph_rad);
544 >                printf("\nback_text alias back_label_mat vellum\n");
545 >                printf("\nback_label_mat polygon back_label\n0\n0\n12");
546 >                printf("\n\t%f %f 0\n\t%f %f 0\n\t%f %f 0\n\t%f %f 0\n",
547 >                                -.25*sph_rad - sph_xoffset, -1.3*sph_rad,
548 >                                -.25*sph_rad - sph_xoffset, (-1.4-1.6*.35/4.-.1)*sph_rad,
549 >                                .25*sph_rad - sph_xoffset, (-1.4-1.6*.35/4.-.1)*sph_rad,
550 >                                .25*sph_rad - sph_xoffset, -1.3*sph_rad );
551 >        }
552 > }
553 >
554 > /* Put out falsecolor scale and name label */
555 > static void
556 > put_scale(void)
557 > {
558 >        const double    max_log10 = log10(overall_max);
559 >        const double    leg_width = 2.*.75*(fabs(sph_xoffset) - sph_rad);
560 >        const double    leg_height = 2.*sph_rad;
561 >        const int       text_lines = 6;
562 >        const int       text_digits = 8;
563 >        char            fmt[16];
564 >        int             i;
565 >
566 >        printf("\n# BSDF legend with falsecolor scale\n");
567 >        printf("\nvoid colorfunc lscale\n10 sca_red(Py) sca_grn(Py) sca_blu(Py)");
568 >        printf("\n\tbsdf2rad.cal -s %f -t 0 %f 0\n0\n0\n", leg_height, -.5*leg_height);
569 >        sprintf(fmt, "%%.%df", text_digits-3);
570 >        for (i = 0; i < text_lines; i++) {
571 >                char    vbuf[16];
572 >                sprintf(vbuf, fmt, pow(10., (i+.5)/text_lines*(max_log10-min_log10)+min_log10));
573 >                printf("\nlscale brighttext lscale\n");
574 >                printf("3 helvet.fnt . %s\n0\n12\n", vbuf);
575 >                printf("\t%f %f 0\n", -.45*leg_width, ((i+.9)/text_lines-.5)*leg_height);
576 >                printf("\t%f 0 0\n", .8*leg_width/strlen(vbuf));
577 >                printf("\t0 %f 0\n", -.9/text_lines*leg_height);
578 >                printf("\t.01 1 -.1\n");
579 >        }
580 >        printf("\nlscale alias legend_mat vellum\n");
581 >        printf("\nlegend_mat polygon legend\n0\n0\n12");
582 >        printf("\n\t%f %f 0\n\t%f %f 0\n\t%f %f 0\n\t%f %f 0\n",
583 >                        -.5*leg_width, .5*leg_height,
584 >                        -.5*leg_width, -.5*leg_height,
585 >                        .5*leg_width, -.5*leg_height,
586 >                        .5*leg_width, .5*leg_height);
587 >        printf("\nvoid brighttext BSDFtitle\n3 helvet.fnt . BSDF\n0\n12\n");
588 >        printf("\t%f %f 0\n", -.25*leg_width, .7*leg_height);
589 >        printf("\t%f 0 0\n", .4/4.*leg_width);
590 >        printf("\t0 %f 0\n", -.1*leg_height);
591 >        printf("\t.01 1 -.1\n");
592 >        printf("\nBSDFtitle alias title_mat vellum\n");
593 >        printf("\ntitle_mat polygon title\n0\n0\n12");
594 >        printf("\n\t%f %f 0\n\t%f %f 0\n\t%f %f 0\n\t%f %f 0\n",
595 >                        -.3*leg_width, .75*leg_height,
596 >                        -.3*leg_width, .55*leg_height,
597 >                        .3*leg_width, .55*leg_height,
598 >                        .3*leg_width, .75*leg_height);
599 >        if (!bsdf_name[0])
600 >                return;
601 >        printf("\nvoid brighttext BSDFname\n3 helvet.fnt . \"%s\"\n0\n12\n", bsdf_name);
602 >        printf("\t%f %f 0\n", -.95*leg_width, -.6*leg_height);
603 >        printf("\t%f 0 0\n", 1.8/strlen(bsdf_name)*leg_width);
604 >        printf("\t0 %f 0\n", -.1*leg_height);
605 >        printf("\t.01 1 -.1\n");
606 >        printf("\nBSDFname alias name_mat vellum\n");
607 >        printf("\nname_mat polygon name\n0\n0\n12");
608 >        printf("\n\t%f %f 0\n\t%f %f 0\n\t%f %f 0\n\t%f %f 0\n",
609 >                        -leg_width, -.55*leg_height,
610 >                        -leg_width, -.75*leg_height,
611 >                        leg_width, -.75*leg_height,
612 >                        leg_width, -.55*leg_height);
613 > }
614 >
615 > /* Convert MGF to Radiance in output */
616 > static void
617 > convert_mgf(const char *mgfdata)
618 > {
619 >        int     len = strlen(mgfdata);
620 >        char    mgfn[128];
621 >        char    radfn[128];
622 >        char    cmdbuf[256];
623 >        float   xmin, xmax, ymin, ymax, zmin, zmax;
624 >        double  max_dim;
625 >        int     fd;
626 >        FILE    *fp;
627 >
628 >        if (!len) return;
629 >        strcpy(mgfn, tfile_name("geom", ".mgf", 0));
630 >        fd = open(mgfn, O_WRONLY|O_CREAT, 0666);
631 >        if (fd < 0 || write(fd, mgfdata, len) != len) {
632 >                fprintf(stderr, "%s: cannot write file '%s'\n",
633 >                                progname, mgfn);
634 >                return;
635 >        }
636 >        close(fd);
637 >        strcpy(radfn, tfile_name("geom", ".rad", 0));
638 >        sprintf(cmdbuf, "mgf2rad %s > %s", mgfn, radfn);
639 >        if (!run_cmd(cmdbuf))
640 >                return;
641 >        sprintf(cmdbuf, "getbbox -w -h %s", radfn);
642 >        if ((fp = popen(cmdbuf, "r")) == NULL ||
643 >                        fscanf(fp, "%f %f %f %f %f %f",
644 >                                &xmin, &xmax, &ymin, &ymax, &zmin, &zmax) != 6
645 >                        || pclose(fp) < 0) {
646 >                fprintf(stderr, "%s: error reading from command: %s\n",
647 >                                progname, cmdbuf);
648 >                return;
649 >        }
650 >        max_dim = ymax - ymin;
651 >        if (xmax - xmin > max_dim)
652 >                max_dim = xmax - xmin;
653 >        if (front_comp) {
654 >                printf("\n# BSDF system geometry (front view)\n");
655 >                sprintf(cmdbuf, "xform -t %f %f %f -s %f -t %f %f 0 %s",
656 >                                -.5*(xmin+xmax), -.5*(ymin+ymax), -zmax,
657 >                                1.5*sph_rad/max_dim,
658 >                                sph_xoffset, -2.5*sph_rad,
659 >                                radfn);
660 >                if (!run_cmd(cmdbuf))
661 >                        return;
662 >        }
663 >        if (back_comp) {
664 >                printf("\n# BSDF system geometry (back view)\n");
665 >                sprintf(cmdbuf, "xform -t %f %f %f -s %f -ry 180 -t %f %f 0 %s",
666 >                                -.5*(xmin+xmax), -.5*(ymin+ymax), -zmin,
667 >                                1.5*sph_rad/max_dim,
668 >                                -sph_xoffset, -2.5*sph_rad,
669 >                                radfn);
670 >                if (!run_cmd(cmdbuf))
671 >                        return;
672 >        }
673 > }
674 >
675 > /* Check RBF input header line & get minimum BSDF value */
676 > static int
677 > rbf_headline(char *s, void *p)
678 > {
679 >        char    fmt[MAXFMTLEN];
680 >
681 >        if (formatval(fmt, s)) {
682 >                if (strcmp(fmt, BSDFREP_FMT))
683 >                        return(-1);
684 >                return(0);
685 >        }
686 >        if (!strncmp(s, "IO_SIDES=", 9)) {
687 >                sscanf(s+9, "%d %d", &input_orient, &output_orient);
688 >                if (input_orient == output_orient) {
689 >                        if (input_orient > 0)
690 >                                front_comp |= SDsampR;
691 >                        else
692 >                                back_comp |= SDsampR;
693 >                } else if (input_orient > 0)
694 >                        front_comp |= SDsampT;
695 >                else
696 >                        back_comp |= SDsampT;
697 >                return(0);
698 >        }
699 >        if (!strncmp(s, "BSDFMIN=", 8)) {
700 >                sscanf(s+8, "%lf", &bsdf_min);
701 >                if (bsdf_min < overall_min)
702 >                        overall_min = bsdf_min;
703 >                return(0);
704 >        }
705 >        return(0);
706 > }
707 >
708 > /* Produce a Radiance model plotting the given BSDF representation */
709   int
710   main(int argc, char *argv[])
711   {
37        int     showPeaks = 0;
38        int     doTrans = 0;
712          int     inpXML = -1;
713 <        RBFNODE *rbf = NULL;
41 <        FILE    *fp;
42 <        char    buf[128];
713 >        double  myLim[2];
714          SDData  myBSDF;
715 <        double  bsdf, min_log;
45 <        FVECT   idir, odir;
46 <        int     i, j, n;
715 >        int     a, n;
716                                                  /* check arguments */
717          progname = argv[0];
718 <        if (argc > 1 && !strcmp(argv[1], "-p")) {
719 <                ++showPeaks;
720 <                ++argv; --argc;
718 >        a = 1;
719 >        myLim[0] = -1; myLim[1] = -2;           /* specified BSDF range? */
720 >        if (argc > a+3 && argv[a][0] == '-' && argv[a][1] == 'r') {
721 >                myLim[0] = atof(argv[++a]);
722 >                myLim[1] = atof(argv[++a]);
723 >                ++a;
724          }
725 <        if (argc > 1 && !strcmp(argv[1], "-t")) {
726 <                ++doTrans;
55 <                ++argv; --argc;
56 <        }
57 <        if (argc >= 4 && (n = strlen(argv[1])-4) > 0) {
58 <                if (!strcasecmp(argv[1]+n, ".xml"))
725 >        if (argc > a && (n = strlen(argv[a])-4) > 0) {
726 >                if (!strcasecmp(argv[a]+n, ".xml"))
727                          inpXML = 1;
728 <                else if (!strcasecmp(argv[1]+n, ".sir"))
728 >                else if (!strcasecmp(argv[a]+n, ".sir"))
729                          inpXML = 0;
730          }
731 <        if (inpXML < 0) {
732 <                fprintf(stderr, "Usage: %s [-p] bsdf.sir theta1 phi1 .. > output.rad\n", progname);
733 <                fprintf(stderr, "   Or: %s [-t] bsdf.xml theta1 phi1 .. > output.rad\n", progname);
731 >        if (inpXML < 0 || inpXML & (argc > a+1)) {
732 >                fprintf(stderr, "Usage: %s [-r min max] bsdf.xml > output.rad\n", progname);
733 >                fprintf(stderr, "   Or: %s [-r min max] hemi1.sir hemi2.sir .. > output.rad\n", progname);
734                  return(1);
735          }
736 <                                                /* load input */
736 >        fputs("# ", stdout);                    /* copy our command */
737 >        printargs(argc, argv, stdout);
738 >                                                /* evaluate BSDF */
739          if (inpXML) {
740 <                SDclearBSDF(&myBSDF, argv[1]);
741 <                if (SDreportError(SDloadFile(&myBSDF, argv[1]), stderr))
740 >                SDclearBSDF(&myBSDF, argv[a]);
741 >                if (SDreportError(SDloadFile(&myBSDF, argv[a]), stderr))
742                          return(1);
743 <                bsdf_min = 1./M_PI;
744 <                if (myBSDF.rf != NULL && myBSDF.rLambFront.cieY < bsdf_min*M_PI)
745 <                        bsdf_min = myBSDF.rLambFront.cieY/M_PI;
746 <                if (myBSDF.rb != NULL && myBSDF.rLambBack.cieY < bsdf_min*M_PI)
747 <                        bsdf_min = myBSDF.rLambBack.cieY/M_PI;
748 <                if ((myBSDF.tf != NULL) | (myBSDF.tb != NULL) &&
749 <                                myBSDF.tLamb.cieY < bsdf_min*M_PI)
80 <                        bsdf_min = myBSDF.tLamb.cieY/M_PI;
81 <                if (doTrans && (myBSDF.tf == NULL) & (myBSDF.tb == NULL)) {
82 <                        fprintf(stderr, "%s: no transmitted component in '%s'\n",
83 <                                        progname, argv[1]);
743 >                if (myBSDF.rf != NULL) front_comp |= SDsampR;
744 >                if (myBSDF.tf != NULL) front_comp |= SDsampT;
745 >                if (myBSDF.rb != NULL) back_comp |= SDsampR;
746 >                if (myBSDF.tb != NULL) back_comp |= SDsampT;
747 >                if (!front_comp & !back_comp) {
748 >                        fprintf(stderr, "%s: nothing to plot in '%s'\n",
749 >                                        progname, argv[a]);
750                          return(1);
751                  }
752 <        } else {
753 <                fp = fopen(argv[1], "rb");
754 <                if (fp == NULL) {
755 <                        fprintf(stderr, "%s: cannot open BSDF interpolant '%s'\n",
756 <                                        progname, argv[1]);
757 <                        return(1);
752 >                if (myLim[0] >= 0)
753 >                        overall_min = myLim[0];
754 >                else {
755 >                        if (front_comp & SDsampR && myBSDF.rLambFront.cieY < overall_min*PI)
756 >                                overall_min = myBSDF.rLambFront.cieY/PI;
757 >                        if (back_comp & SDsampR && myBSDF.rLambBack.cieY < overall_min*PI)
758 >                                overall_min = myBSDF.rLambBack.cieY/PI;
759 >                        if ((front_comp|back_comp) & SDsampT &&
760 >                                        myBSDF.tLamb.cieY < overall_min*PI)
761 >                                overall_min = myBSDF.tLamb.cieY/PI;
762                  }
763 <                if (!load_bsdf_rep(fp))
763 >                set_minlog();
764 >                if (!build_wBSDF(&myBSDF))
765                          return(1);
766 <                fclose(fp);
767 <        }
768 < #ifdef DEBUG
769 <        fprintf(stderr, "Minimum BSDF set to %.4f\n", bsdf_min);
770 < #endif
771 <        min_log = log(bsdf_min*.5 + 1e-5);
772 <                                                /* output BSDF rep. */
773 <        for (n = 0; (n < 6) & (2*n+3 < argc); n++) {
774 <                double  theta = (M_PI/180.)*atof(argv[2*n+2]);
775 <                double  phi = (M_PI/180.)*atof(argv[2*n+3]);
105 <                if (theta < -FTINY) {
106 <                        fprintf(stderr, "%s: theta values must be positive\n",
107 <                                        progname);
766 >                if (myBSDF.matn[0])
767 >                        strcpy(bsdf_name, myBSDF.matn);
768 >                else
769 >                        strcpy(bsdf_name, myBSDF.name);
770 >                strcpy(bsdf_manuf, myBSDF.makr);
771 >                put_matBSDF(argv[a]);
772 >        } else {
773 >                FILE    *fp[4];
774 >                if (argc > a+4) {
775 >                        fprintf(stderr, "%s: more than 4 hemispheres!\n", progname);
776                          return(1);
777                  }
778 <                if (inpXML) {
779 <                        input_orient = (theta <= M_PI/2.) ? 1 : -1;
780 <                        output_orient = doTrans ? -input_orient : input_orient;
778 >                for (n = a; n < argc; n++) {
779 >                        fp[n-a] = fopen(argv[n], "rb");
780 >                        if (fp[n-a] == NULL) {
781 >                                fprintf(stderr, "%s: cannot open BSDF interpolant '%s'\n",
782 >                                                progname, argv[n]);
783 >                                return(1);
784 >                        }
785 >                        if (getheader(fp[n-a], rbf_headline, NULL) < 0) {
786 >                                fprintf(stderr, "%s: bad BSDF interpolant '%s'\n",
787 >                                                progname, argv[n]);
788 >                                return(1);
789 >                        }
790                  }
791 <                idir[2] = sin(theta);
792 <                idir[0] = idir[2] * cos(phi);
793 <                idir[1] = idir[2] * sin(phi);
794 <                idir[2] = input_orient * sqrt(1. - idir[2]*idir[2]);
795 < #ifdef DEBUG
796 <                fprintf(stderr, "Computing BSDF for incident direction (%.1f,%.1f)\n",
797 <                                get_theta180(idir), get_phi360(idir));
798 < #endif
122 <                if (!inpXML)
123 <                        rbf = advect_rbf(idir, 15000);
124 < #ifdef DEBUG
125 <                if (inpXML)
126 <                        fprintf(stderr, "Hemispherical %s: %.3f\n",
127 <                                (output_orient > 0 ^ input_orient > 0 ?
128 <                                        "transmission" : "reflection"),
129 <                                SDdirectHemi(idir, SDsampSp|SDsampDf |
130 <                                        (output_orient > 0 ^ input_orient > 0 ?
131 <                                                 SDsampT : SDsampR), &myBSDF));
132 <                else if (rbf == NULL)
133 <                        fputs("Empty RBF\n", stderr);
134 <                else
135 <                        fprintf(stderr, "Hemispherical %s: %.3f\n",
136 <                                (output_orient > 0 ^ input_orient > 0 ?
137 <                                        "transmission" : "reflection"),
138 <                                rbf->vtotal);
139 < #endif
140 <                printf("void trans tmat\n0\n0\n7 %f %f %f .04 .04 .9 1\n",
141 <                                colarr[n][0], colarr[n][1], colarr[n][2]);
142 <                if (showPeaks && rbf != NULL) {
143 <                        printf("void plastic pmat\n0\n0\n5 %f %f %f .04 .08\n",
144 <                                1.-colarr[n][0], 1.-colarr[n][1], 1.-colarr[n][2]);
145 <                        for (i = 0; i < rbf->nrbf; i++) {
146 <                                ovec_from_pos(odir, rbf->rbfa[i].gx, rbf->rbfa[i].gy);
147 <                                bsdf = eval_rbfrep(rbf, odir);
148 <                                bsdf = log(bsdf + 1e-5) - min_log;
149 <                                printf("pmat sphere p%d\n0\n0\n4 %f %f %f %f\n",
150 <                                        i+1, odir[0]*bsdf, odir[1]*bsdf, odir[2]*bsdf,
151 <                                                .007*bsdf);
791 >                if (myLim[0] >= 0)
792 >                        overall_min = myLim[0];
793 >                set_minlog();
794 >                for (n = a; n < argc; n++) {
795 >                        if (fseek(fp[n-a], 0L, SEEK_SET) < 0) {
796 >                                fprintf(stderr, "%s: cannot seek on '%s'\n",
797 >                                                progname, argv[n]);
798 >                                return(1);
799                          }
800 +                        if (!load_bsdf_rep(fp[n-a]))
801 +                                return(1);
802 +                        fclose(fp[n-a]);
803 +                        if (!build_wRBF())
804 +                                return(1);
805                  }
806 <                fflush(stdout);
155 <                sprintf(buf, "gensurf tmat bsdf%d - - - %d %d %s", n+1,
156 <                                                GRIDRES-1, GRIDRES-1, validf);
157 <                fp = popen(buf, "w");
158 <                if (fp == NULL) {
159 <                        fprintf(stderr, "%s: cannot open '| %s'\n", progname, buf);
160 <                        return(1);
161 <                }
162 <                for (i = 0; i < GRIDRES; i++)
163 <                    for (j = 0; j < GRIDRES; j++) {
164 <                        ovec_from_pos(odir, i, j);
165 <                        if (inpXML) {
166 <                                SDValue sval;
167 <                                if (SDreportError(SDevalBSDF(&sval, odir,
168 <                                                        idir, &myBSDF), stderr))
169 <                                        return(1);
170 <                                bsdf = sval.cieY;
171 <                        } else
172 <                                bsdf = eval_rbfrep(rbf, odir);
173 <                        bsdf = log(bsdf + 1e-5) - min_log;
174 <                        fprintf(fp, "%.8e %.8e %.8e\n",
175 <                                        odir[0]*bsdf, odir[1]*bsdf, odir[2]*bsdf);
176 <                    }
177 <                if (rbf != NULL)
178 <                        free(rbf);
179 <                if (pclose(fp))
180 <                        return(1);
806 >                put_matBSDF(NULL);
807          }
808 +        if (myLim[1] > myLim[0])        /* override maximum BSDF? */
809 +                overall_max = myLim[1];
810 +        put_source();                   /* before hemispheres & labels */
811 +        put_hemispheres();
812 +        put_scale();
813 +        if (inpXML && myBSDF.mgf)
814 +                convert_mgf(myBSDF.mgf);
815 +        if (!put_BSDFs())               /* most of the output happens here */
816 +                return(1);
817 +        cleanup_tmp();
818          return(0);
819   }

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines