ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/cv/bsdf2rado.c
Revision: 2.5
Committed: Tue Dec 7 23:55:02 2021 UTC (2 years, 4 months ago) by greg
Content type: text/plain
Branch: MAIN
CVS Tags: rad5R4, HEAD
Changes since 2.4: +3 -3 lines
Log Message:
fix: Repaired issue with reciprocity and BSDF sampling, thanks to Jon Sargent

File Contents

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