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

Comparing ray/src/cv/bsdf2klems.c (file contents):
Revision 2.9 by greg, Sun Aug 11 14:32:34 2013 UTC vs.
Revision 2.24 by greg, Fri Feb 17 22:31:49 2017 UTC

# Line 14 | Line 14 | static const char RCSid[] = "$Id$";
14   #include <math.h>
15   #include "random.h"
16   #include "platform.h"
17 + #include "paths.h"
18 + #include "rtio.h"
19   #include "calcomp.h"
20   #include "bsdfrep.h"
21   #include "bsdf_m.h"
22 +                                /* tristimulus components */
23 + enum {CIE_X, CIE_Y, CIE_Z};
24                                  /* assumed maximum # Klems patches */
25   #define MAXPATCHES      145
26                                  /* global argv[0] */
27   char                    *progname;
28                                  /* selected basis function name */
29 < static const char       *kbasis = "LBNL/Klems Full";
29 > static const char       klems_full[] = "LBNL/Klems Full";
30 > static const char       klems_half[] = "LBNL/Klems Half";
31 > static const char       klems_quarter[] = "LBNL/Klems Quarter";
32 > static const char       *kbasis = klems_full;
33                                  /* number of BSDF samples per patch */
34   static int              npsamps = 256;
35 +                                /* limit on number of RBF lobes */
36 + static int              lobe_lim = 15000;
37 +                                /* progress bar length */
38 + static int              do_prog = 79;
39  
40 < /* Return angle basis corresponding to the given name */
41 < ANGLE_BASIS *
42 < get_basis(const char *bn)
32 < {
33 <        int     n = nabases;
34 <        
35 <        while (n-- > 0)
36 <                if (!strcasecmp(bn, abase_list[n].name))
37 <                        return &abase_list[n];
38 <        return NULL;
39 < }
40 > #define MAXCARG         512     /* wrapBSDF command */
41 > static char             *wrapBSDF[MAXCARG] = {"wrapBSDF", "-W", "-UU"};
42 > static int              wbsdfac = 3;
43  
44 < /* Output XML header to stdout */
44 > /* Add argument to wrapBSDF, allocating space if !isstatic */
45   static void
46 < xml_header(int ac, char *av[])
46 > add_wbsdf(const char *arg, int isstatic)
47   {
48 <        puts("<?xml version=\"1.0\" encoding=\"UTF-8\"?>");
49 <        puts("<WindowElement xmlns=\"http://windows.lbl.gov\" xmlns:xsi=\"http://www.w3.org/2001/XMLSchema-instance\" xsi:schemaLocation=\"http://windows.lbl.gov/BSDF-v1.4.xsd\">");
50 <        fputs("<!-- File produced by:", stdout);
51 <        while (ac-- > 0) {
52 <                fputc(' ', stdout);
53 <                fputs(*av++, stdout);
48 >        if (arg == NULL)
49 >                return;
50 >        if (wbsdfac >= MAXCARG-1) {
51 >                fputs(progname, stderr);
52 >                fputs(": too many command arguments to wrapBSDF\n", stderr);
53 >                exit(1);
54          }
55 <        puts(" -->");
55 >        if (!*arg)
56 >                arg = "";
57 >        else if (!isstatic)
58 >                arg = savqstr((char *)arg);
59 >
60 >        wrapBSDF[wbsdfac++] = (char *)arg;
61   }
62  
63 < /* Output XML prologue to stdout */
63 > /* Start new progress bar */
64 > #define prog_start(s)   if (do_prog) fprintf(stderr, "%s: %s...\n", progname, s); else
65 >
66 > /* Draw progress bar of the appropriate length */
67   static void
68 < xml_prologue(const SDData *sd)
68 > prog_show(double frac)
69   {
70 <        const char      *matn = (sd && sd->matn[0]) ? sd->matn : "Name";
71 <        const char      *makr = (sd && sd->makr[0]) ? sd->makr : "Manufacturer";
72 <        ANGLE_BASIS     *abp = get_basis(kbasis);
73 <        int             i;
70 >        static unsigned call_cnt = 0;
71 >        static char     lastc[] = "-\\|/";
72 >        char            pbar[256];
73 >        int             nchars;
74  
75 <        if (abp == NULL) {
76 <                fprintf(stderr, "%s: unknown angle basis '%s'\n", progname, kbasis);
77 <                exit(1);
78 <        }
79 <        puts("<WindowElementType>System</WindowElementType>");
80 <        puts("<FileType>BSDF</FileType>");
81 <        puts("<Optical>");
82 <        puts("<Layer>");
83 <        puts("\t<Material>");
84 <        printf("\t\t<Name>%s</Name>\n", matn);
85 <        printf("\t\t<Manufacturer>%s</Manufacturer>\n", makr);
86 <        if (sd && sd->dim[2] > .001) {
76 <                printf("\t\t<Thickness unit=\"meter\">%.3f</Thickness>\n", sd->dim[2]);
77 <                printf("\t\t<Width unit=\"meter\">%.3f</Width>\n", sd->dim[0]);
78 <                printf("\t\t<Height unit=\"meter\">%.3f</Height>\n", sd->dim[1]);
79 <        }
80 <        puts("\t\t<DeviceType>Other</DeviceType>");
81 <        puts("\t</Material>");
82 <        if (sd && sd->mgf != NULL) {
83 <                puts("\t<Geometry format=\"MGF\">");
84 <                puts("\t\t<MGFblock unit=\"meter\">");
85 <                fputs(sd->mgf, stdout);
86 <                puts("</MGFblock>");
87 <                puts("\t</Geometry>");
88 <        }
89 <        puts("\t<DataDefinition>");
90 <        puts("\t\t<IncidentDataStructure>Columns</IncidentDataStructure>");
91 <        puts("\t\t<AngleBasis>");
92 <        printf("\t\t\t<AngleBasisName>%s</AngleBasisName>\n", kbasis);
93 <        for (i = 0; abp->lat[i].nphis; i++) {
94 <                puts("\t\t\t<AngleBasisBlock>");
95 <                printf("\t\t\t<Theta>%g</Theta>\n", i ?
96 <                                .5*(abp->lat[i].tmin + abp->lat[i+1].tmin) :
97 <                                .0 );
98 <                printf("\t\t\t<nPhis>%d</nPhis>\n", abp->lat[i].nphis);
99 <                puts("\t\t\t<ThetaBounds>");
100 <                printf("\t\t\t\t<LowerTheta>%g</LowerTheta>\n", abp->lat[i].tmin);
101 <                printf("\t\t\t\t<UpperTheta>%g</UpperTheta>\n", abp->lat[i+1].tmin);
102 <                puts("\t\t\t</ThetaBounds>");
103 <                puts("\t\t\t</AngleBasisBlock>");
104 <        }
105 <        puts("\t\t</AngleBasis>");
106 <        puts("\t</DataDefinition>");
75 >        if (do_prog <= 1) return;
76 >        if (do_prog > sizeof(pbar)-2)
77 >                do_prog = sizeof(pbar)-2;
78 >        if (frac < 0) frac = 0;
79 >        else if (frac >= 1) frac = .9999;
80 >        nchars = do_prog*frac;
81 >        pbar[0] = '\r';
82 >        memset(pbar+1, '*', nchars);
83 >        pbar[nchars+1] = lastc[call_cnt++ & 3];
84 >        memset(pbar+2+nchars, '-', do_prog-nchars-1);
85 >        pbar[do_prog+1] = '\0';
86 >        fputs(pbar, stderr);
87   }
88  
89 < /* Output XML data prologue to stdout */
89 > /* Finish progress bar */
90   static void
91 < data_prologue()
91 > prog_done(void)
92   {
93 <        static const char       *bsdf_type[4] = {
114 <                                        "Reflection Front",
115 <                                        "Transmission Front",
116 <                                        "Transmission Back",
117 <                                        "Reflection Back"
118 <                                };
93 >        int     n = do_prog;
94  
95 <        puts("\t<WavelengthData>");
96 <        puts("\t\t<LayerNumber>System</LayerNumber>");
97 <        puts("\t\t<Wavelength unit=\"Integral\">Visible</Wavelength>");
98 <        puts("\t\t<SourceSpectrum>CIE Illuminant D65 1nm.ssp</SourceSpectrum>");
99 <        puts("\t\t<DetectorSpectrum>ASTM E308 1931 Y.dsp</DetectorSpectrum>");
125 <        puts("\t\t<WavelengthDataBlock>");
126 <        printf("\t\t\t<WavelengthDataDirection>%s</WavelengthDataDirection>\n",
127 <                        bsdf_type[(input_orient>0)<<1 | (output_orient>0)]);
128 <        printf("\t\t\t<ColumnAngleBasis>%s</ColumnAngleBasis>\n", kbasis);
129 <        printf("\t\t\t<RowAngleBasis>%s</RowAngleBasis>\n", kbasis);
130 <        puts("\t\t\t<ScatteringDataType>BTDF</ScatteringDataType>");
131 <        puts("\t\t\t<ScatteringData>");
95 >        if (n <= 1) return;
96 >        fputc('\r', stderr);
97 >        while (n--)
98 >                fputc(' ', stderr);
99 >        fputc('\r', stderr);
100   }
101  
102 < /* Output XML data epilogue to stdout */
103 < static void
104 < data_epilogue(void)
102 > /* Return angle basis corresponding to the given name */
103 > static ANGLE_BASIS *
104 > get_basis(const char *bn)
105   {
106 <        puts("\t\t\t</ScatteringData>");
107 <        puts("\t\t</WavelengthDataBlock>");
108 <        puts("\t</WavelengthData>");
106 >        int     n = nabases;
107 >        
108 >        while (n-- > 0)
109 >                if (!strcasecmp(bn, abase_list[n].name))
110 >                        return &abase_list[n];
111 >        return NULL;
112   }
113  
114 < /* Output XML epilogue to stdout */
115 < static void
116 < xml_epilogue(void)
114 > /* Copy geometry string to file for wrapBSDF */
115 > static char *
116 > save_geom(const char *mgf)
117   {
118 <        puts("</Layer>");
119 <        puts("</Optical>");
120 <        puts("</WindowElement>");
118 >        char    *tfname = mktemp(savqstr(TEMPLATE));
119 >        int     fd = open(tfname, O_CREAT|O_WRONLY, 0600);
120 >
121 >        if (fd < 0)
122 >                return(NULL);
123 >        write(fd, mgf, strlen(mgf));
124 >        close(fd);
125 >        add_wbsdf("-g", 1);
126 >        add_wbsdf(tfname, 1);
127 >        return(tfname);
128   }
129  
130 + /* Open XYZ component file for output and add appropriate arguments */
131 + static FILE *
132 + open_component_file(int c)
133 + {
134 +        static const char       sname[3][6] = {"CIE-X", "CIE-Y", "CIE-Z"};
135 +        static const char       cname[4][4] = {"-rf", "-tf", "-tb", "-rb"};
136 +        char                    *tfname = mktemp(savqstr(TEMPLATE));
137 +        FILE                    *fp = fopen(tfname, "w");
138 +
139 +        if (fp == NULL) {
140 +                fprintf(stderr, "%s: cannot open '%s' for writing\n",
141 +                                progname, tfname);
142 +                exit(1);
143 +        }
144 +        add_wbsdf("-s", 1); add_wbsdf(sname[c], 1);
145 +        add_wbsdf(cname[(input_orient>0)<<1 | (output_orient>0)], 1);
146 +        add_wbsdf(tfname, 1);
147 +        return(fp);
148 + }
149 +
150   /* Load and resample XML BSDF description using Klems basis */
151   static void
152   eval_bsdf(const char *fname)
153   {
154          ANGLE_BASIS     *abp = get_basis(kbasis);
155 +        FILE            *cfp[3];
156          SDData          bsd;
157          SDError         ec;
158          FVECT           vin, vout;
159 <        SDValue         sv;
160 <        double          sum;
159 >        SDValue         sdv;
160 >        double          sum, xsum, ysum;
161          int             i, j, n;
162  
163 +        initurand(npsamps);
164          SDclearBSDF(&bsd, fname);               /* load BSDF file */
165          if ((ec = SDloadFile(&bsd, fname)) != SDEnone)
166                  goto err;
167 <        xml_prologue(&bsd);                     /* pass geometry */
167 >        if (bsd.mgf != NULL)                    /* save geometry */
168 >                save_geom(bsd.mgf);
169 >        if (bsd.matn[0])                        /* save identifier(s) */
170 >                strcpy(bsdf_name, bsd.matn);
171 >        if (bsd.makr[0])
172 >                strcpy(bsdf_manuf, bsd.makr);
173 >        if (bsd.dim[2] > 0) {                   /* save dimension(s) */
174 >                char    buf[64];
175 >                if ((bsd.dim[0] > 0) & (bsd.dim[1] > 0))
176 >                        sprintf(buf, "w=%g;h=%g;t=%g",
177 >                                        bsd.dim[0], bsd.dim[1], bsd.dim[2]);
178 >                else
179 >                        sprintf(buf, "t=%g", bsd.dim[2]);
180 >                add_wbsdf("-f", 1);
181 >                add_wbsdf(buf, 0);
182 >        }
183                                                  /* front reflection */
184          if (bsd.rf != NULL || bsd.rLambFront.cieY > .002) {
185              input_orient = 1; output_orient = 1;
186 <            data_prologue();
186 >            cfp[CIE_Y] = open_component_file(CIE_Y);
187 >            if (bsd.rf != NULL && bsd.rf->comp[0].cspec[2].flags) {
188 >                rbf_colorimetry = RBCtristimulus;
189 >                cfp[CIE_X] = open_component_file(CIE_X);
190 >                cfp[CIE_Z] = open_component_file(CIE_Z);
191 >            } else
192 >                rbf_colorimetry = RBCphotopic;
193              for (j = 0; j < abp->nangles; j++) {
194                  for (i = 0; i < abp->nangles; i++) {
195                      sum = 0;                    /* average over patches */
196 +                    xsum = ysum = 0;
197                      for (n = npsamps; n-- > 0; ) {
198                          fo_getvec(vout, j+(n+frandom())/npsamps, abp);
199                          fi_getvec(vin, i+urand(n), abp);
200 <                        ec = SDevalBSDF(&sv, vout, vin, &bsd);
200 >                        ec = SDevalBSDF(&sdv, vout, vin, &bsd);
201                          if (ec != SDEnone)
202                                  goto err;
203 <                        sum += sv.cieY;
203 >                        sum += sdv.cieY;
204 >                        if (rbf_colorimetry == RBCtristimulus) {
205 >                                xsum += sdv.cieY * sdv.spec.cx;
206 >                                ysum += sdv.cieY * sdv.spec.cy;
207 >                        }
208                      }
209 <                    printf("\t%.3e\n", sum/npsamps);
209 >                    fprintf(cfp[CIE_Y], "\t%.3e\n", sum/npsamps);
210 >                    if (rbf_colorimetry == RBCtristimulus) {
211 >                        fprintf(cfp[CIE_X], "\t%.3e\n", xsum*sum/(npsamps*ysum));
212 >                        fprintf(cfp[CIE_Z], "\t%.3e\n",
213 >                                (sum - xsum - ysum)*sum/(npsamps*ysum));
214 >                    }
215                  }
216 <                putchar('\n');                  /* extra space between rows */
216 >                fputc('\n', cfp[CIE_Y]);        /* extra space between rows */
217 >                if (rbf_colorimetry == RBCtristimulus) {
218 >                        fputc('\n', cfp[CIE_X]);
219 >                        fputc('\n', cfp[CIE_Z]);
220 >                }
221              }
222 <            data_epilogue();
222 >            if (fclose(cfp[CIE_Y])) {
223 >                fprintf(stderr, "%s: error writing Y output\n", progname);
224 >                exit(1);
225 >            }
226 >            if (rbf_colorimetry == RBCtristimulus &&
227 >                        (fclose(cfp[CIE_X]) || fclose(cfp[CIE_Z]))) {
228 >                fprintf(stderr, "%s: error writing X/Z output\n", progname);
229 >                exit(1);
230 >            }
231          }
232                                                  /* back reflection */
233          if (bsd.rb != NULL || bsd.rLambBack.cieY > .002) {
234              input_orient = -1; output_orient = -1;
235 <            data_prologue();
235 >            cfp[CIE_Y] = open_component_file(CIE_Y);
236 >            if (bsd.rb != NULL && bsd.rb->comp[0].cspec[2].flags) {
237 >                rbf_colorimetry = RBCtristimulus;
238 >                cfp[CIE_X] = open_component_file(CIE_X);
239 >                cfp[CIE_Z] = open_component_file(CIE_Z);
240 >            } else
241 >                rbf_colorimetry = RBCphotopic;
242              for (j = 0; j < abp->nangles; j++) {
243                  for (i = 0; i < abp->nangles; i++) {
244                      sum = 0;                    /* average over patches */
245 +                    xsum = ysum = 0;
246                      for (n = npsamps; n-- > 0; ) {
247                          bo_getvec(vout, j+(n+frandom())/npsamps, abp);
248                          bi_getvec(vin, i+urand(n), abp);
249 <                        ec = SDevalBSDF(&sv, vout, vin, &bsd);
249 >                        ec = SDevalBSDF(&sdv, vout, vin, &bsd);
250                          if (ec != SDEnone)
251                                  goto err;
252 <                        sum += sv.cieY;
252 >                        sum += sdv.cieY;
253 >                        if (rbf_colorimetry == RBCtristimulus) {
254 >                                xsum += sdv.cieY * sdv.spec.cx;
255 >                                ysum += sdv.cieY * sdv.spec.cy;
256 >                        }
257                      }
258 <                    printf("\t%.3e\n", sum/npsamps);
258 >                    fprintf(cfp[CIE_Y], "\t%.3e\n", sum/npsamps);
259 >                    if (rbf_colorimetry == RBCtristimulus) {
260 >                        fprintf(cfp[CIE_X], "\t%.3e\n", xsum*sum/(npsamps*ysum));
261 >                        fprintf(cfp[CIE_Z], "\t%.3e\n",
262 >                                (sum - xsum - ysum)*sum/(npsamps*ysum));
263 >                    }
264                  }
265 <                putchar('\n');                  /* extra space between rows */
265 >                if (rbf_colorimetry == RBCtristimulus) {
266 >                        fputc('\n', cfp[CIE_X]);
267 >                        fputc('\n', cfp[CIE_Z]);
268 >                }
269              }
270 <            data_epilogue();
270 >            if (fclose(cfp[CIE_Y])) {
271 >                fprintf(stderr, "%s: error writing Y output\n", progname);
272 >                exit(1);
273 >            }
274 >            if (rbf_colorimetry == RBCtristimulus &&
275 >                        (fclose(cfp[CIE_X]) || fclose(cfp[CIE_Z]))) {
276 >                fprintf(stderr, "%s: error writing X/Z output\n", progname);
277 >                exit(1);
278 >            }
279          }
280                                                  /* front transmission */
281          if (bsd.tf != NULL || bsd.tLamb.cieY > .002) {
282              input_orient = 1; output_orient = -1;
283 <            data_prologue();
283 >            cfp[CIE_Y] = open_component_file(CIE_Y);
284 >            if (bsd.tf != NULL && bsd.tf->comp[0].cspec[2].flags) {
285 >                rbf_colorimetry = RBCtristimulus;
286 >                cfp[CIE_X] = open_component_file(CIE_X);
287 >                cfp[CIE_Z] = open_component_file(CIE_Z);
288 >            } else
289 >                rbf_colorimetry = RBCphotopic;
290              for (j = 0; j < abp->nangles; j++) {
291                  for (i = 0; i < abp->nangles; i++) {
292                      sum = 0;                    /* average over patches */
293 +                    xsum = ysum = 0;
294                      for (n = npsamps; n-- > 0; ) {
295                          bo_getvec(vout, j+(n+frandom())/npsamps, abp);
296                          fi_getvec(vin, i+urand(n), abp);
297 <                        ec = SDevalBSDF(&sv, vout, vin, &bsd);
297 >                        ec = SDevalBSDF(&sdv, vout, vin, &bsd);
298                          if (ec != SDEnone)
299                                  goto err;
300 <                        sum += sv.cieY;
300 >                        sum += sdv.cieY;
301 >                        if (rbf_colorimetry == RBCtristimulus) {
302 >                                xsum += sdv.cieY * sdv.spec.cx;
303 >                                ysum += sdv.cieY * sdv.spec.cy;
304 >                        }
305                      }
306 <                    printf("\t%.3e\n", sum/npsamps);
306 >                    fprintf(cfp[CIE_Y], "\t%.3e\n", sum/npsamps);
307 >                    if (rbf_colorimetry == RBCtristimulus) {
308 >                        fprintf(cfp[CIE_X], "\t%.3e\n", xsum*sum/(npsamps*ysum));
309 >                        fprintf(cfp[CIE_Z], "\t%.3e\n",
310 >                                (sum - xsum - ysum)*sum/(npsamps*ysum));
311 >                    }
312                  }
313 <                putchar('\n');                  /* extra space between rows */
313 >                if (rbf_colorimetry == RBCtristimulus) {
314 >                        fputc('\n', cfp[CIE_X]);
315 >                        fputc('\n', cfp[CIE_Z]);
316 >                }
317              }
318 <            data_epilogue();
318 >            if (fclose(cfp[CIE_Y])) {
319 >                fprintf(stderr, "%s: error writing Y output\n", progname);
320 >                exit(1);
321 >            }
322 >            if (rbf_colorimetry == RBCtristimulus &&
323 >                        (fclose(cfp[CIE_X]) || fclose(cfp[CIE_Z]))) {
324 >                fprintf(stderr, "%s: error writing X/Z output\n", progname);
325 >                exit(1);
326 >            }
327          }
328                                                  /* back transmission */
329          if ((bsd.tb != NULL) | (bsd.tf != NULL)) {
330              input_orient = -1; output_orient = 1;
331 <            data_prologue();
331 >            cfp[CIE_Y] = open_component_file(CIE_Y);
332 >            if (bsd.tb != NULL)
333 >                rbf_colorimetry = bsd.tb->comp[0].cspec[2].flags
334 >                                        ? RBCtristimulus : RBCphotopic ;
335 >            if (rbf_colorimetry == RBCtristimulus) {
336 >                cfp[CIE_X] = open_component_file(CIE_X);
337 >                cfp[CIE_Z] = open_component_file(CIE_Z);
338 >            }
339              for (j = 0; j < abp->nangles; j++) {
340                  for (i = 0; i < abp->nangles; i++) {
341                      sum = 0;            /* average over patches */
342 +                    xsum = ysum = 0;
343                      for (n = npsamps; n-- > 0; ) {
344                          fo_getvec(vout, j+(n+frandom())/npsamps, abp);
345                          bi_getvec(vin, i+urand(n), abp);
346 <                        ec = SDevalBSDF(&sv, vout, vin, &bsd);
346 >                        ec = SDevalBSDF(&sdv, vout, vin, &bsd);
347                          if (ec != SDEnone)
348                                  goto err;
349 <                        sum += sv.cieY;
349 >                        sum += sdv.cieY;
350 >                        if (rbf_colorimetry == RBCtristimulus) {
351 >                                xsum += sdv.cieY * sdv.spec.cx;
352 >                                ysum += sdv.cieY * sdv.spec.cy;
353 >                        }
354                      }
355 <                    printf("\t%.3e\n", sum/npsamps);
355 >                    fprintf(cfp[CIE_Y], "\t%.3e\n", sum/npsamps);
356 >                    if (rbf_colorimetry == RBCtristimulus) {
357 >                        fprintf(cfp[CIE_X], "\t%.3e\n", xsum*sum/(npsamps*ysum));
358 >                        fprintf(cfp[CIE_Z], "\t%.3e\n",
359 >                                (sum - xsum - ysum)*sum/(npsamps*ysum));
360 >                    }
361                  }
362 <                putchar('\n');                  /* extra space between rows */
362 >                if (rbf_colorimetry == RBCtristimulus) {
363 >                        fputc('\n', cfp[CIE_X]);
364 >                        fputc('\n', cfp[CIE_Z]);
365 >                }
366              }
367 <            data_epilogue();
367 >            if (fclose(cfp[CIE_Y])) {
368 >                fprintf(stderr, "%s: error writing Y output\n", progname);
369 >                exit(1);
370 >            }
371 >            if (rbf_colorimetry == RBCtristimulus &&
372 >                        (fclose(cfp[CIE_X]) || fclose(cfp[CIE_Z]))) {
373 >                fprintf(stderr, "%s: error writing X/Z output\n", progname);
374 >                exit(1);
375 >            }
376          }
377          SDfreeBSDF(&bsd);                       /* all done */
378          return;
# Line 262 | Line 387 | eval_function(char *funame)
387   {
388          ANGLE_BASIS     *abp = get_basis(kbasis);
389          int             assignD = (fundefined(funame) < 6);
390 +        FILE            *ofp = open_component_file(CIE_Y);
391          double          iovec[6];
392          double          sum;
393          int             i, j, n;
394  
395          initurand(npsamps);
270        data_prologue();                        /* begin output */
396          for (j = 0; j < abp->nangles; j++) {    /* run through directions */
397              for (i = 0; i < abp->nangles; i++) {
398                  sum = 0;
# Line 290 | Line 415 | eval_function(char *funame)
415                      }
416                      sum += funvalue(funame, 6, iovec);
417                  }
418 <                printf("\t%.3e\n", sum/npsamps);
418 >                fprintf(ofp, "\t%.3e\n", sum/npsamps);
419              }
420 <            putchar('\n');
420 >            fputc('\n', ofp);
421 >            prog_show((j+1.)/abp->nangles);
422          }
423 <        data_epilogue();                        /* finish output */
423 >        prog_done();
424 >        if (fclose(ofp)) {
425 >                fprintf(stderr, "%s: error writing Y output\n", progname);
426 >                exit(1);
427 >        }
428   }
429  
430   /* Interpolate and output a radial basis function BSDF representation */
# Line 302 | Line 432 | static void
432   eval_rbf(void)
433   {
434          ANGLE_BASIS     *abp = get_basis(kbasis);
435 +        float           (*XZarr)[2] = NULL;
436          float           bsdfarr[MAXPATCHES*MAXPATCHES];
437 +        FILE            *cfp[3];
438          FVECT           vin, vout;
439 <        RBFNODE         *rbf;
308 <        double          sum;
439 >        double          sum, xsum, ysum;
440          int             i, j, n;
441                                                  /* sanity check */
442          if (abp->nangles > MAXPATCHES) {
443                  fprintf(stderr, "%s: too many patches!\n", progname);
444                  exit(1);
445          }
446 <        data_prologue();                        /* begin output */
446 >        if (rbf_colorimetry == RBCtristimulus)
447 >                XZarr = (float (*)[2])malloc(sizeof(float)*2*abp->nangles*abp->nangles);
448          for (i = 0; i < abp->nangles; i++) {
449 +            RBFNODE     *rbf;
450              if (input_orient > 0)               /* use incident patch center */
451                  fi_getvec(vin, i+.5*(i>0), abp);
452              else
453                  bi_getvec(vin, i+.5*(i>0), abp);
454  
455 <            rbf = advect_rbf(vin);              /* compute radial basis func */
455 >            rbf = advect_rbf(vin, lobe_lim);    /* compute radial basis func */
456  
457              for (j = 0; j < abp->nangles; j++) {
458                  sum = 0;                        /* sample over exiting patch */
459 +                xsum = ysum = 0;
460                  for (n = npsamps; n--; ) {
461 +                    SDValue     sdv;
462                      if (output_orient > 0)
463                          fo_getvec(vout, j+(n+frandom())/npsamps, abp);
464                      else
465                          bo_getvec(vout, j+(n+frandom())/npsamps, abp);
466  
467 <                    sum += eval_rbfrep(rbf, vout) / vout[2];
467 >                    eval_rbfcol(&sdv, rbf, vout);
468 >                    sum += sdv.cieY;
469 >                    if (rbf_colorimetry == RBCtristimulus) {
470 >                        xsum += sdv.cieY * sdv.spec.cx;
471 >                        ysum += sdv.cieY * sdv.spec.cy;
472 >                    }
473                  }
474 <                bsdfarr[j*abp->nangles + i] = sum*output_orient/npsamps;
474 >                n = j*abp->nangles + i;
475 >                bsdfarr[n] = sum / npsamps;
476 >                if (rbf_colorimetry == RBCtristimulus) {
477 >                    XZarr[n][0] = xsum*sum/(npsamps*ysum);
478 >                    XZarr[n][1] = (sum - xsum - ysum)*sum/(npsamps*ysum);
479 >                }
480              }
481              if (rbf != NULL)
482                  free(rbf);
483 +            prog_show((i+1.)/abp->nangles);
484          }
485 <        n = 0;                                  /* write out our matrix */
485 >                                                /* write out our matrix */
486 >        cfp[CIE_Y] = open_component_file(CIE_Y);
487 >        n = 0;
488          for (j = 0; j < abp->nangles; j++) {
489 <            for (i = 0; i < abp->nangles; i++)
490 <                printf("\t%.3e\n", bsdfarr[n++]);
491 <            putchar('\n');
489 >            for (i = 0; i < abp->nangles; i++, n++)
490 >                fprintf(cfp[CIE_Y], "\t%.3e\n", bsdfarr[n]);
491 >            fputc('\n', cfp[CIE_Y]);
492          }
493 <        data_epilogue();                        /* finish output */
493 >        prog_done();
494 >        if (fclose(cfp[CIE_Y])) {
495 >                fprintf(stderr, "%s: error writing Y output\n", progname);
496 >                exit(1);
497 >        }
498 >        if (XZarr == NULL)                      /* no color? */
499 >                return;
500 >        cfp[CIE_X] = open_component_file(CIE_X);
501 >        cfp[CIE_Z] = open_component_file(CIE_Z);
502 >        n = 0;
503 >        for (j = 0; j < abp->nangles; j++) {
504 >            for (i = 0; i < abp->nangles; i++, n++) {
505 >                fprintf(cfp[CIE_X], "\t%.3e\n", XZarr[n][0]);
506 >                fprintf(cfp[CIE_Z], "\t%.3e\n", XZarr[n][1]);
507 >            }
508 >            fputc('\n', cfp[CIE_X]);
509 >            fputc('\n', cfp[CIE_Z]);
510 >        }
511 >        free(XZarr);
512 >        if (fclose(cfp[CIE_X]) || fclose(cfp[CIE_Z])) {
513 >                fprintf(stderr, "%s: error writing X/Z output\n", progname);
514 >                exit(1);
515 >        }
516   }
517  
518 + #if defined(_WIN32) || defined(_WIN64)
519 + /* Execute wrapBSDF command (may never return) */
520 + static int
521 + wrap_up(void)
522 + {
523 +        char    cmd[8192];
524 +
525 +        if (bsdf_manuf[0]) {
526 +                add_wbsdf("-f", 1);
527 +                strcpy(cmd, "m=");
528 +                strcpy(cmd+2, bsdf_manuf);
529 +                add_wbsdf(cmd, 0);
530 +        }
531 +        if (bsdf_name[0]) {
532 +                add_wbsdf("-f", 1);
533 +                strcpy(cmd, "n=");
534 +                strcpy(cmd+2, bsdf_name);
535 +                add_wbsdf(cmd, 0);
536 +        }
537 +        if (!convert_commandline(cmd, sizeof(cmd), wrapBSDF)) {
538 +                fputs(progname, stderr);
539 +                fputs(": command line too long in wrap_up()\n", stderr);
540 +                return(1);
541 +        }
542 +        return(system(cmd));
543 + }
544 + #else
545 + /* Execute wrapBSDF command (may never return) */
546 + static int
547 + wrap_up(void)
548 + {
549 +        char    buf[256];
550 +        char    *compath = getpath((char *)wrapBSDF[0], getenv("PATH"), X_OK);
551 +
552 +        if (compath == NULL) {
553 +                fprintf(stderr, "%s: cannot locate %s\n", progname, wrapBSDF[0]);
554 +                return(1);
555 +        }
556 +        if (bsdf_manuf[0]) {
557 +                add_wbsdf("-f", 1);
558 +                strcpy(buf, "m=");
559 +                strcpy(buf+2, bsdf_manuf);
560 +                add_wbsdf(buf, 0);
561 +        }
562 +        if (bsdf_name[0]) {
563 +                add_wbsdf("-f", 1);
564 +                strcpy(buf, "n=");
565 +                strcpy(buf+2, bsdf_name);
566 +                add_wbsdf(buf, 0);
567 +        }
568 +        execv(compath, wrapBSDF);       /* successful call never returns */
569 +        perror(compath);
570 +        return(1);
571 + }
572 + #endif
573 +
574   /* Read in BSDF and interpolate as Klems matrix representation */
575   int
576   main(int argc, char *argv[])
577   {
578          int     dofwd = 0, dobwd = 1;
579 +        char    buf[2048];
580          char    *cp;
581          int     i, na;
582  
# Line 371 | Line 598 | main(int argc, char *argv[])
598                          break;
599                  case 'f':
600                          if (!argv[i][2]) {
601 <                                fcompile(argv[++i]);
602 <                                single_plane_incident = 0;
601 >                                if (strchr(argv[++i], '=') != NULL) {
602 >                                        add_wbsdf("-f", 1);
603 >                                        add_wbsdf(argv[i], 1);
604 >                                } else {
605 >                                        fcompile(argv[i]);
606 >                                        single_plane_incident = 0;
607 >                                }
608                          } else
609                                  dofwd = (argv[i][0] == '+');
610                          break;
# Line 380 | Line 612 | main(int argc, char *argv[])
612                          dobwd = (argv[i][0] == '+');
613                          break;
614                  case 'h':
615 <                        kbasis = "LBNL/Klems Half";
615 >                        kbasis = klems_half;
616 >                        add_wbsdf("-a", 1);
617 >                        add_wbsdf("kh", 1);
618                          break;
619                  case 'q':
620 <                        kbasis = "LBNL/Klems Quarter";
620 >                        kbasis = klems_quarter;
621 >                        add_wbsdf("-a", 1);
622 >                        add_wbsdf("kq", 1);
623                          break;
624 +                case 'l':
625 +                        lobe_lim = atoi(argv[++i]);
626 +                        break;
627 +                case 'p':
628 +                        do_prog = atoi(argv[i]+2);
629 +                        break;
630 +                case 'C':
631 +                        add_wbsdf(argv[i], 1);
632 +                        add_wbsdf(argv[++i], 1);
633 +                        break;
634                  default:
635                          goto userr;
636                  }
637 +        if (kbasis == klems_full) {             /* default (full) basis? */
638 +                add_wbsdf("-a", 1);
639 +                add_wbsdf("kf", 1);
640 +        }
641 +        strcpy(buf, "File produced by: ");
642 +        if (convert_commandline(buf+18, sizeof(buf)-18, argv) != NULL) {
643 +                add_wbsdf("-C", 1); add_wbsdf(buf, 0);
644 +        }
645          if (single_plane_incident >= 0) {       /* function-based BSDF? */
646 <                if (i != argc-1 || fundefined(argv[i]) != 6) {
646 >                if (i != argc-1 || fundefined(argv[i]) < 3) {
647                          fprintf(stderr,
648          "%s: need single function with 6 arguments: bsdf(ix,iy,iz,ox,oy,oz)\n",
649                                          progname);
650 <                        fprintf(stderr, "\tor 3 arguments using Dx,Dy,Dz: bsdf(ix,iy,iz)\n",
397 <                                        progname);
650 >                        fprintf(stderr, "\tor 3 arguments using Dx,Dy,Dz: bsdf(ix,iy,iz)\n");
651                          goto userr;
652                  }
653                  ++eclock;
401                xml_header(argc, argv);                 /* start XML output */
402                xml_prologue(NULL);
654                  if (dofwd) {
655                          input_orient = -1;
656                          output_orient = -1;
657 <                        eval_function(argv[i]);         /* outside reflectance */
657 >                        prog_start("Evaluating outside reflectance");
658 >                        eval_function(argv[i]);
659                          output_orient = 1;
660 <                        eval_function(argv[i]);         /* outside -> inside */
660 >                        prog_start("Evaluating outside->inside transmission");
661 >                        eval_function(argv[i]);
662                  }
663                  if (dobwd) {
664                          input_orient = 1;
665                          output_orient = 1;
666 <                        eval_function(argv[i]);         /* inside reflectance */
666 >                        prog_start("Evaluating inside reflectance");
667 >                        eval_function(argv[i]);
668                          output_orient = -1;
669 <                        eval_function(argv[i]);         /* inside -> outside */
669 >                        prog_start("Evaluating inside->outside transmission");
670 >                        eval_function(argv[i]);
671                  }
672 <                xml_epilogue();                 /* finish XML output & exit */
418 <                return(0);
672 >                return(wrap_up());
673          }
674                                                  /* XML input? */
675          if (i == argc-1 && (cp = argv[i]+strlen(argv[i])-4) > argv[i] &&
676                                  !strcasecmp(cp, ".xml")) {
423                xml_header(argc, argv);         /* start XML output */
677                  eval_bsdf(argv[i]);             /* load & resample BSDF */
678 <                xml_epilogue();                 /* finish XML output & exit */
426 <                return(0);
678 >                return(wrap_up());
679          }
680          if (i < argc) {                         /* open input files if given */
681                  int     nbsdf = 0;
682                  for ( ; i < argc; i++) {        /* interpolate each component */
683 +                        char    pbuf[256];
684                          FILE    *fpin = fopen(argv[i], "rb");
685                          if (fpin == NULL) {
686                                  fprintf(stderr, "%s: cannot open BSDF interpolant '%s'\n",
# Line 437 | Line 690 | main(int argc, char *argv[])
690                          if (!load_bsdf_rep(fpin))
691                                  return(1);
692                          fclose(fpin);
693 <                        if (!nbsdf++) {         /* start XML on first dist. */
694 <                                xml_header(argc, argv);
442 <                                xml_prologue(NULL);
443 <                        }
693 >                        sprintf(pbuf, "Interpolating component '%s'", argv[i]);
694 >                        prog_start(pbuf);
695                          eval_rbf();
696                  }
697 <                xml_epilogue();                 /* finish XML output & exit */
447 <                return(0);
697 >                return(wrap_up());
698          }
699          SET_FILE_BINARY(stdin);                 /* load from stdin */
700          if (!load_bsdf_rep(stdin))
701                  return(1);
702 <        xml_header(argc, argv);                 /* start XML output */
453 <        xml_prologue(NULL);
702 >        prog_start("Interpolating from standard input");
703          eval_rbf();                             /* resample dist. */
704 <        xml_epilogue();                         /* finish XML output & exit */
456 <        return(0);
704 >        return(wrap_up());
705   userr:
706          fprintf(stderr,
707 <        "Usage: %s [-n spp][-h|-q][bsdf.sir ..] > bsdf.xml\n", progname);
707 >        "Usage: %s [-n spp][-h|-q][-l maxlobes] [bsdf.sir ..] > bsdf.xml\n", progname);
708          fprintf(stderr,
709          "   or: %s [-n spp][-h|-q] bsdf_in.xml > bsdf_out.xml\n", progname);
710          fprintf(stderr,

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines