ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/cv/bsdf2rad.c
Revision: 2.10
Committed: Sun Mar 30 17:23:28 2014 UTC (10 years, 1 month ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 2.9: +5 -3 lines
Log Message:
Fixed handling of zero BSDF function evaluations

File Contents

# Content
1 #ifndef lint
2 static const char RCSid[] = "$Id: bsdf2rad.c,v 2.9 2014/03/24 06:07:46 greg Exp $";
3 #endif
4 /*
5 * Plot 3-D BSDF output based on scattering interpolant or XML representation
6 */
7
8 #define _USE_MATH_DEFINES
9 #include <stdio.h>
10 #include <string.h>
11 #include <stdlib.h>
12 #include <math.h>
13 #include "rtprocess.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 };
24
25 char validf[] = "-e \"valid(s,t)=X`SYS(s,t)^2+Y`SYS(s,t)^2+Z`SYS(s,t)^2-1e-7\"";
26
27 char *progname;
28
29 /* Produce a Radiance model plotting the indicated incident direction(s) */
30 int
31 main(int argc, char *argv[])
32 {
33 int showPeaks = 0;
34 int doTrans = 0;
35 int inpXML = -1;
36 RBFNODE *rbf = NULL;
37 FILE *fp;
38 char buf[128];
39 SDData myBSDF;
40 double bsdf, min_log;
41 FVECT idir, odir;
42 int i, j, n;
43 /* check arguments */
44 progname = argv[0];
45 if (argc > 1 && !strcmp(argv[1], "-p")) {
46 ++showPeaks;
47 ++argv; --argc;
48 }
49 if (argc > 1 && !strcmp(argv[1], "-t")) {
50 ++doTrans;
51 ++argv; --argc;
52 }
53 if (argc >= 4 && (n = strlen(argv[1])-4) > 0) {
54 if (!strcasecmp(argv[1]+n, ".xml"))
55 inpXML = 1;
56 else if (!strcasecmp(argv[1]+n, ".sir"))
57 inpXML = 0;
58 }
59 if (inpXML < 0) {
60 fprintf(stderr, "Usage: %s [-p] bsdf.sir theta1 phi1 .. > output.rad\n", progname);
61 fprintf(stderr, " Or: %s [-t] bsdf.xml theta1 phi1 .. > output.rad\n", progname);
62 return(1);
63 }
64 /* load input */
65 if (inpXML) {
66 SDclearBSDF(&myBSDF, argv[1]);
67 if (SDreportError(SDloadFile(&myBSDF, argv[1]), stderr))
68 return(1);
69 bsdf_min = 1./M_PI;
70 if (myBSDF.rf != NULL && myBSDF.rLambFront.cieY < bsdf_min*M_PI)
71 bsdf_min = myBSDF.rLambFront.cieY/M_PI;
72 if (myBSDF.rb != NULL && myBSDF.rLambBack.cieY < bsdf_min*M_PI)
73 bsdf_min = myBSDF.rLambBack.cieY/M_PI;
74 if ((myBSDF.tf != NULL) | (myBSDF.tb != NULL) &&
75 myBSDF.tLamb.cieY < bsdf_min*M_PI)
76 bsdf_min = myBSDF.tLamb.cieY/M_PI;
77 if (doTrans && (myBSDF.tf == NULL) & (myBSDF.tb == NULL)) {
78 fprintf(stderr, "%s: no transmitted component in '%s'\n",
79 progname, argv[1]);
80 return(1);
81 }
82 } else {
83 fp = fopen(argv[1], "rb");
84 if (fp == NULL) {
85 fprintf(stderr, "%s: cannot open BSDF interpolant '%s'\n",
86 progname, argv[1]);
87 return(1);
88 }
89 if (!load_bsdf_rep(fp))
90 return(1);
91 fclose(fp);
92 }
93 #ifdef DEBUG
94 fprintf(stderr, "Minimum BSDF set to %.4f\n", bsdf_min);
95 #endif
96 min_log = log(bsdf_min*.5 + 1e-5);
97 /* output BSDF rep. */
98 for (n = 0; (n < 6) & (2*n+3 < argc); n++) {
99 double theta = (M_PI/180.)*atof(argv[2*n+2]);
100 double phi = (M_PI/180.)*atof(argv[2*n+3]);
101 if (theta < -FTINY) {
102 fprintf(stderr, "%s: theta values must be positive\n",
103 progname);
104 return(1);
105 }
106 if (inpXML) {
107 input_orient = (theta <= M_PI/2.) ? 1 : -1;
108 output_orient = doTrans ? -input_orient : input_orient;
109 }
110 idir[2] = sin(theta);
111 idir[0] = idir[2] * cos(phi);
112 idir[1] = idir[2] * sin(phi);
113 idir[2] = input_orient * sqrt(1. - idir[2]*idir[2]);
114 #ifdef DEBUG
115 fprintf(stderr, "Computing BSDF for incident direction (%.1f,%.1f)\n",
116 get_theta180(idir), get_phi360(idir));
117 #endif
118 if (!inpXML)
119 rbf = advect_rbf(idir, 15000);
120 #ifdef DEBUG
121 if (inpXML)
122 fprintf(stderr, "Hemispherical %s: %.3f\n",
123 (output_orient > 0 ? "reflection" : "transmission"),
124 SDdirectHemi(idir, SDsampSp|SDsampDf |
125 (output_orient > 0 ?
126 SDsampR : SDsampT), &myBSDF));
127 else if (rbf == NULL)
128 fputs("Empty RBF\n", stderr);
129 else
130 fprintf(stderr, "Hemispherical %s: %.3f\n",
131 (output_orient > 0 ? "reflection" : "transmission"),
132 rbf->vtotal);
133 #endif
134 printf("void trans tmat\n0\n0\n7 %f %f %f .04 .04 .9 1\n",
135 colarr[n][0], colarr[n][1], colarr[n][2]);
136 if (showPeaks && rbf != NULL) {
137 printf("void plastic pmat\n0\n0\n5 %f %f %f .04 .08\n",
138 1.-colarr[n][0], 1.-colarr[n][1], 1.-colarr[n][2]);
139 for (i = 0; i < rbf->nrbf; i++) {
140 ovec_from_pos(odir, rbf->rbfa[i].gx, rbf->rbfa[i].gy);
141 bsdf = eval_rbfrep(rbf, odir) / (output_orient*odir[2]);
142 bsdf = log(bsdf + 1e-5) - min_log;
143 printf("pmat sphere p%d\n0\n0\n4 %f %f %f %f\n",
144 i+1, odir[0]*bsdf, odir[1]*bsdf, odir[2]*bsdf,
145 .007*bsdf);
146 }
147 }
148 fflush(stdout);
149 sprintf(buf, "gensurf tmat bsdf%d - - - %d %d %s", n+1,
150 GRIDRES-1, GRIDRES-1, validf);
151 fp = popen(buf, "w");
152 if (fp == NULL) {
153 fprintf(stderr, "%s: cannot open '| %s'\n", progname, buf);
154 return(1);
155 }
156 for (i = 0; i < GRIDRES; i++)
157 for (j = 0; j < GRIDRES; j++) {
158 ovec_from_pos(odir, i, j);
159 if (inpXML) {
160 SDValue sval;
161 if (SDreportError(SDevalBSDF(&sval, odir,
162 idir, &myBSDF), stderr))
163 return(1);
164 bsdf = sval.cieY;
165 } else
166 bsdf = eval_rbfrep(rbf, odir) /
167 (output_orient*odir[2]);
168 bsdf = log(bsdf + 1e-5) - min_log;
169 fprintf(fp, "%.8e %.8e %.8e\n",
170 odir[0]*bsdf, odir[1]*bsdf, odir[2]*bsdf);
171 }
172 if (rbf != NULL)
173 free(rbf);
174 if (pclose(fp))
175 return(1);
176 }
177 return(0);
178 }