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

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines