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

Comparing ray/src/cv/bsdf2ttree.c (file contents):
Revision 2.13 by greg, Sat Mar 23 02:21:14 2013 UTC vs.
Revision 2.32 by greg, Tue Feb 2 18:02:32 2016 UTC

# Line 11 | Line 11 | static const char RCSid[] = "$Id$";
11   #include <stdio.h>
12   #include <stdlib.h>
13   #include <math.h>
14 + #include "random.h"
15   #include "platform.h"
16 + #include "rtprocess.h"
17   #include "calcomp.h"
18   #include "bsdfrep.h"
19                                  /* global argv[0] */
20   char                    *progname;
21                                  /* percentage to cull (<0 to turn off) */
22 < double                  pctcull = 90.;
22 > static double           pctcull = 90.;
23                                  /* sampling order */
24 < int                     samp_order = 6;
24 > static int              samp_order = 6;
25 >                                /* super-sampling threshold */
26 > const double            ssamp_thresh = 0.35;
27 >                                /* number of super-samples */
28 > #ifndef NSSAMP
29 > #define NSSAMP          100
30 > #endif
31 >                                /* limit on number of RBF lobes */
32 > static int              lobe_lim = 15000;
33 >                                /* progress bar length */
34 > static int              do_prog = 79;
35  
24 /* Output XML prologue to stdout */
25 static void
26 xml_prologue(int ac, char *av[])
27 {
28        puts("<?xml version=\"1.0\" encoding=\"UTF-8\"?>");
29        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\">");
30        fputs("<!-- File produced by:", stdout);
31        while (ac-- > 0) {
32                fputc(' ', stdout);
33                fputs(*av++, stdout);
34        }
35        puts(" -->");
36        puts("<WindowElementType>System</WindowElementType>");
37        puts("<FileType>BSDF</FileType>");
38        puts("<Optical>");
39        puts("<Layer>");
40        puts("\t<Material>");
41        puts("\t\t<Name>Name</Name>");
42        puts("\t\t<Manufacturer>Manufacturer</Manufacturer>");
43        puts("\t\t<DeviceType>Other</DeviceType>");
44        puts("\t</Material>");
45        puts("\t<DataDefinition>");
46        printf("\t\t<IncidentDataStructure>TensorTree%c</IncidentDataStructure>\n",
47                        single_plane_incident ? '3' : '4');
48        puts("\t</DataDefinition>");
49 }
36  
37 < /* Output XML data prologue to stdout */
37 >
38 > /* Start new progress bar */
39 > #define prog_start(s)   if (do_prog) fprintf(stderr, "%s: %s...\n", progname, s); else
40 >
41 > /* Draw progress bar of the appropriate length */
42   static void
43 < data_prologue()
43 > prog_show(double frac)
44   {
45 <        static const char       *bsdf_type[4] = {
46 <                                        "Reflection Front",
47 <                                        "Transmission Front",
48 <                                        "Transmission Back",
59 <                                        "Reflection Back"
60 <                                };
45 >        static unsigned call_cnt = 0;
46 >        static char     lastc[] = "-\\|/";
47 >        char            pbar[256];
48 >        int             nchars;
49  
50 <        puts("\t<WavelengthData>");
51 <        puts("\t\t<LayerNumber>System</LayerNumber>");
52 <        puts("\t\t<Wavelength unit=\"Integral\">Visible</Wavelength>");
53 <        puts("\t\t<SourceSpectrum>CIE Illuminant D65 1nm.ssp</SourceSpectrum>");
54 <        puts("\t\t<DetectorSpectrum>ASTM E308 1931 Y.dsp</DetectorSpectrum>");
55 <        puts("\t\t<WavelengthDataBlock>");
56 <        printf("\t\t\t<WavelengthDataDirection>%s</WavelengthDataDirection>\n",
57 <                        bsdf_type[(input_orient>0)<<1 | (output_orient>0)]);
58 <        puts("\t\t\t<AngleBasis>LBNL/Shirley-Chiu</AngleBasis>");
59 <        puts("\t\t\t<ScatteringDataType>BTDF</ScatteringDataType>");
60 <        puts("\t\t\t<ScatteringData>");
50 >        if (do_prog <= 1) return;
51 >        if (do_prog > sizeof(pbar)-2)
52 >                do_prog = sizeof(pbar)-2;
53 >        if (frac < 0) frac = 0;
54 >        else if (frac >= 1) frac = .9999;
55 >        nchars = do_prog*frac;
56 >        pbar[0] = '\r';
57 >        memset(pbar+1, '*', nchars);
58 >        pbar[nchars+1] = lastc[call_cnt++ & 3];
59 >        memset(pbar+2+nchars, '-', do_prog-nchars-1);
60 >        pbar[do_prog+1] = '\0';
61 >        fputs(pbar, stderr);
62   }
63  
64 < /* Output XML data epilogue to stdout */
64 > /* Finish progress bar */
65   static void
66 < data_epilogue(void)
66 > prog_done(void)
67   {
68 <        puts("\t\t\t</ScatteringData>");
69 <        puts("\t\t</WavelengthDataBlock>");
70 <        puts("\t</WavelengthData>");
68 >        int     n = do_prog;
69 >
70 >        if (n <= 1) return;
71 >        fputc('\r', stderr);
72 >        while (n--)
73 >                fputc(' ', stderr);
74 >        fputc('\r', stderr);
75   }
76  
77 < /* Output XML epilogue to stdout */
78 < static void
79 < xml_epilogue(void)
77 > /* Compute absolute relative difference */
78 > static double
79 > abs_diff(double v1, double v0)
80   {
81 <        puts("</Layer>");
82 <        puts("</Optical>");
83 <        puts("</WindowElement>");
81 >        if ((v0 < 0) | (v1 < 0))
82 >                return(.0);
83 >        v1 = (v1-v0)*2./(v0+v1+.0001);
84 >        if (v1 < 0)
85 >                return(-v1);
86 >        return(v1);
87   }
88  
89   /* Interpolate and output isotropic BSDF data */
# Line 96 | Line 92 | eval_isotropic(char *funame)
92   {
93          const int       sqres = 1<<samp_order;
94          FILE            *ofp = NULL;
95 +        int             assignD = 0;
96          char            cmd[128];
97          int             ix, ox, oy;
98          double          iovec[6];
99          float           bsdf;
100 < #if DEBUG
104 <        fprintf(stderr, "Writing isotropic order %d ", samp_order);
105 <        if (pctcull >= 0) fprintf(stderr, "data with %.1f%% culling\n", pctcull);
106 <        else fputs("raw data\n", stderr);
107 < #endif
100 >
101          data_prologue();                        /* begin output */
102          if (pctcull >= 0) {
103 <                sprintf(cmd, "rttree_reduce -h -a -ff -r 3 -t %f -g %d",
103 >                sprintf(cmd, "rttree_reduce -a -h -ff -r 3 -t %f -g %d",
104                                  pctcull, samp_order);
105                  fflush(stdout);
106                  ofp = popen(cmd, "w");
# Line 117 | Line 110 | eval_isotropic(char *funame)
110                          exit(1);
111                  }
112                  SET_FILE_BINARY(ofp);
113 + #ifdef getc_unlocked                            /* avoid lock/unlock overhead */
114 +                flockfile(ofp);
115 + #endif
116          } else
117                  fputs("{\n", stdout);
118 +                                                /* need to assign Dx, Dy, Dz? */
119 +        if (funame != NULL)
120 +                assignD = (fundefined(funame) < 6);
121                                                  /* run through directions */
122          for (ix = 0; ix < sqres/2; ix++) {
123                  RBFNODE *rbf = NULL;
124 <                iovec[0] = (ix+.5)/sqres - 1.;
124 >                iovec[0] = 2.*(ix+.5)/sqres - 1.;
125                  iovec[1] = .0;
126                  iovec[2] = input_orient * sqrt(1. - iovec[0]*iovec[0]);
127                  if (funame == NULL)
128 <                        rbf = advect_rbf(iovec);
129 <                for (ox = 0; ox < sqres; ox++)
128 >                        rbf = advect_rbf(iovec, lobe_lim);
129 >                for (ox = 0; ox < sqres; ox++) {
130 >                    float       last_bsdf = -1;
131                      for (oy = 0; oy < sqres; oy++) {
132                          SDsquare2disk(iovec+3, (ox+.5)/sqres, (oy+.5)/sqres);
133                          iovec[5] = output_orient *
134                                  sqrt(1. - iovec[3]*iovec[3] - iovec[4]*iovec[4]);
135                          if (funame == NULL)
136 <                                bsdf = eval_rbfrep(rbf, iovec+3) *
137 <                                                output_orient/iovec[5];
138 <                        else
139 <                                bsdf = funvalue(funame, 6, iovec);
136 >                            bsdf = eval_rbfrep(rbf, iovec+3);
137 >                        else {
138 >                            if (assignD) {
139 >                                varset("Dx", '=', -iovec[3]);
140 >                                varset("Dy", '=', -iovec[4]);
141 >                                varset("Dz", '=', -iovec[5]);
142 >                                ++eclock;
143 >                            }
144 >                            bsdf = funvalue(funame, 6, iovec);
145 > #if (NSSAMP > 0)
146 >                            if (abs_diff(bsdf, last_bsdf) > ssamp_thresh) {
147 >                                int     ssi;
148 >                                double  ssa[3], ssvec[6], sum = 0;
149 >                                                /* super-sample voxel */
150 >                                for (ssi = NSSAMP; ssi--; ) {
151 >                                    SDmultiSamp(ssa, 3, (ssi+frandom()) *
152 >                                                        (1./NSSAMP));
153 >                                    ssvec[0] = 2.*(ix+ssa[0])/sqres - 1.;
154 >                                    ssvec[1] = .0;
155 >                                    ssvec[2] = input_orient *
156 >                                                sqrt(1. - ssvec[0]*ssvec[0]);
157 >                                    SDsquare2disk(ssvec+3, (ox+ssa[1])/sqres,
158 >                                                (oy+ssa[2])/sqres);
159 >                                    ssvec[5] = output_orient *
160 >                                                sqrt(1. - ssvec[3]*ssvec[3] -
161 >                                                        ssvec[4]*ssvec[4]);
162 >                                    if (assignD) {
163 >                                        varset("Dx", '=', -ssvec[3]);
164 >                                        varset("Dy", '=', -ssvec[4]);
165 >                                        varset("Dz", '=', -ssvec[5]);
166 >                                        ++eclock;
167 >                                    }
168 >                                    sum += funvalue(funame, 6, ssvec);
169 >                                }
170 >                                bsdf = sum/NSSAMP;
171 >                            }
172 > #endif
173 >                        }
174                          if (pctcull >= 0)
175                                  fwrite(&bsdf, sizeof(bsdf), 1, ofp);
176                          else
177                                  printf("\t%.3e\n", bsdf);
178 +                        last_bsdf = bsdf;
179                      }
180 +                }
181                  if (rbf != NULL)
182                          free(rbf);
183 +                prog_show((ix+1.)*(2./sqres));
184          }
185          if (pctcull >= 0) {                     /* finish output */
186                  if (pclose(ofp)) {
# Line 157 | Line 194 | eval_isotropic(char *funame)
194                  fputs("}\n", stdout);
195          }
196          data_epilogue();
197 +        prog_done();
198   }
199  
200   /* Interpolate and output anisotropic BSDF data */
# Line 165 | Line 203 | eval_anisotropic(char *funame)
203   {
204          const int       sqres = 1<<samp_order;
205          FILE            *ofp = NULL;
206 +        int             assignD = 0;
207          char            cmd[128];
208          int             ix, iy, ox, oy;
209          double          iovec[6];
210          float           bsdf;
211 < #if DEBUG
173 <        fprintf(stderr, "Writing anisotropic order %d ", samp_order);
174 <        if (pctcull >= 0) fprintf(stderr, "data with %.1f%% culling\n", pctcull);
175 <        else fputs("raw data\n", stderr);
176 < #endif
211 >
212          data_prologue();                        /* begin output */
213          if (pctcull >= 0) {
214 <                sprintf(cmd, "rttree_reduce -h -a -ff -r 4 -t %f -g %d",
214 >                sprintf(cmd, "rttree_reduce%s -h -ff -r 4 -t %f -g %d",
215 >                                (input_orient>0 ^ output_orient>0) ? "" : " -a",
216                                  pctcull, samp_order);
217                  fflush(stdout);
218                  ofp = popen(cmd, "w");
# Line 185 | Line 221 | eval_anisotropic(char *funame)
221                                          progname);
222                          exit(1);
223                  }
224 +                SET_FILE_BINARY(ofp);
225 + #ifdef getc_unlocked                            /* avoid lock/unlock overhead */
226 +                flockfile(ofp);
227 + #endif
228          } else
229                  fputs("{\n", stdout);
230 +                                                /* need to assign Dx, Dy, Dz? */
231 +        if (funame != NULL)
232 +                assignD = (fundefined(funame) < 6);
233                                                  /* run through directions */
234          for (ix = 0; ix < sqres; ix++)
235              for (iy = 0; iy < sqres; iy++) {
236                  RBFNODE *rbf = NULL;            /* Klems reversal */
237 <                SDsquare2disk(iovec, (ix+.5)/sqres, (iy+.5)/sqres);
195 <                iovec[0] = -iovec[0]; iovec[1] = -iovec[1];
237 >                SDsquare2disk(iovec, 1.-(ix+.5)/sqres, 1.-(iy+.5)/sqres);
238                  iovec[2] = input_orient *
239                                  sqrt(1. - iovec[0]*iovec[0] - iovec[1]*iovec[1]);
240                  if (funame == NULL)
241 <                        rbf = advect_rbf(iovec);
242 <                for (ox = 0; ox < sqres; ox++)
241 >                        rbf = advect_rbf(iovec, lobe_lim);
242 >                for (ox = 0; ox < sqres; ox++) {
243 >                    float       last_bsdf = -1;
244                      for (oy = 0; oy < sqres; oy++) {
245                          SDsquare2disk(iovec+3, (ox+.5)/sqres, (oy+.5)/sqres);
246                          iovec[5] = output_orient *
247                                  sqrt(1. - iovec[3]*iovec[3] - iovec[4]*iovec[4]);
248                          if (funame == NULL)
249 <                                bsdf = eval_rbfrep(rbf, iovec+3) *
250 <                                                output_orient/iovec[5];
251 <                        else
252 <                                bsdf = funvalue(funame, 6, iovec);
249 >                            bsdf = eval_rbfrep(rbf, iovec+3);
250 >                        else {
251 >                            if (assignD) {
252 >                                varset("Dx", '=', -iovec[3]);
253 >                                varset("Dy", '=', -iovec[4]);
254 >                                varset("Dz", '=', -iovec[5]);
255 >                                ++eclock;
256 >                            }
257 >                            bsdf = funvalue(funame, 6, iovec);
258 > #if (NSSAMP > 0)
259 >                            if (abs_diff(bsdf, last_bsdf) > ssamp_thresh) {
260 >                                int     ssi;
261 >                                double  ssa[4], ssvec[6], sum = 0;
262 >                                                /* super-sample voxel */
263 >                                for (ssi = NSSAMP; ssi--; ) {
264 >                                    SDmultiSamp(ssa, 4, (ssi+frandom()) *
265 >                                                        (1./NSSAMP));
266 >                                    SDsquare2disk(ssvec, 1.-(ix+ssa[0])/sqres,
267 >                                                1.-(iy+ssa[1])/sqres);
268 >                                    ssvec[2] = input_orient *
269 >                                                sqrt(1. - ssvec[0]*ssvec[0] -
270 >                                                        ssvec[1]*ssvec[1]);
271 >                                    SDsquare2disk(ssvec+3, (ox+ssa[2])/sqres,
272 >                                                (oy+ssa[3])/sqres);
273 >                                    ssvec[5] = output_orient *
274 >                                                sqrt(1. - ssvec[3]*ssvec[3] -
275 >                                                        ssvec[4]*ssvec[4]);
276 >                                    if (assignD) {
277 >                                        varset("Dx", '=', -ssvec[3]);
278 >                                        varset("Dy", '=', -ssvec[4]);
279 >                                        varset("Dz", '=', -ssvec[5]);
280 >                                        ++eclock;
281 >                                    }
282 >                                    sum += funvalue(funame, 6, ssvec);
283 >                                }
284 >                                bsdf = sum/NSSAMP;
285 >                            }
286 > #endif
287 >                        }
288                          if (pctcull >= 0)
289                                  fwrite(&bsdf, sizeof(bsdf), 1, ofp);
290                          else
291                                  printf("\t%.3e\n", bsdf);
292 +                        last_bsdf = bsdf;
293                      }
294 +                }
295                  if (rbf != NULL)
296                          free(rbf);
297 +                prog_show((ix*sqres+iy+1.)/(sqres*sqres));
298              }
299          if (pctcull >= 0) {                     /* finish output */
300                  if (pclose(ofp)) {
# Line 224 | Line 305 | eval_anisotropic(char *funame)
305          } else
306                  fputs("}\n", stdout);
307          data_epilogue();
308 +        prog_done();
309   }
310  
311   /* Read in BSDF and interpolate as tensor tree representation */
# Line 270 | Line 352 | main(int argc, char *argv[])
352                  case 'g':
353                          samp_order = atoi(argv[++i]);
354                          break;
355 +                case 'l':
356 +                        lobe_lim = atoi(argv[++i]);
357 +                        break;
358 +                case 'p':
359 +                        do_prog = atoi(argv[i]+2);
360 +                        break;
361                  default:
362                          goto userr;
363                  }
364          if (single_plane_incident >= 0) {       /* function-based BSDF? */
365                  void    (*evf)(char *s) = single_plane_incident ?
366                                  &eval_isotropic : &eval_anisotropic;
367 <                if (i != argc-1 || fundefined(argv[i]) != 6) {
367 >                if (i != argc-1 || fundefined(argv[i]) < 3) {
368                          fprintf(stderr,
369          "%s: need single function with 6 arguments: bsdf(ix,iy,iz,ox,oy,oz)\n",
370                                          progname);
371 +                        fprintf(stderr, "\tor 3 arguments using Dx,Dy,Dz: bsdf(ix,iy,iz)\n");
372                          goto userr;
373                  }
374 +                ++eclock;
375                  xml_prologue(argc, argv);       /* start XML output */
376                  if (dofwd) {
377                          input_orient = -1;
378                          output_orient = -1;
379 <                        (*evf)(argv[i]);        /* outside reflectance */
379 >                        prog_start("Evaluating outside reflectance");
380 >                        (*evf)(argv[i]);
381                          output_orient = 1;
382 <                        (*evf)(argv[i]);        /* outside -> inside */
382 >                        prog_start("Evaluating outside->inside transmission");
383 >                        (*evf)(argv[i]);
384                  }
385                  if (dobwd) {
386                          input_orient = 1;
387                          output_orient = 1;
388 <                        (*evf)(argv[i]);        /* inside reflectance */
388 >                        prog_start("Evaluating inside reflectance");
389 >                        (*evf)(argv[i]);
390                          output_orient = -1;
391 <                        (*evf)(argv[i]);        /* inside -> outside */
391 >                        prog_start("Evaluating inside->outside transmission");
392 >                        (*evf)(argv[i]);
393                  }
394                  xml_epilogue();                 /* finish XML output & exit */
395                  return(0);
# Line 303 | Line 397 | main(int argc, char *argv[])
397          if (i < argc) {                         /* open input files if given */
398                  int     nbsdf = 0;
399                  for ( ; i < argc; i++) {        /* interpolate each component */
400 +                        char    pbuf[256];
401                          FILE    *fpin = fopen(argv[i], "rb");
402                          if (fpin == NULL) {
403                                  fprintf(stderr, "%s: cannot open BSDF interpolant '%s'\n",
# Line 314 | Line 409 | main(int argc, char *argv[])
409                          fclose(fpin);
410                          if (!nbsdf++)           /* start XML on first dist. */
411                                  xml_prologue(argc, argv);
412 +                        sprintf(pbuf, "Interpolating component '%s'", argv[i]);
413 +                        prog_start(pbuf);
414                          if (single_plane_incident)
415                                  eval_isotropic(NULL);
416                          else
# Line 326 | Line 423 | main(int argc, char *argv[])
423          if (!load_bsdf_rep(stdin))
424                  return(1);
425          xml_prologue(argc, argv);               /* start XML output */
426 +        prog_start("Interpolating from standard input");
427          if (single_plane_incident)              /* resample dist. */
428                  eval_isotropic(NULL);
429          else
# Line 334 | Line 432 | main(int argc, char *argv[])
432          return(0);
433   userr:
434          fprintf(stderr,
435 <        "Usage: %s [-g Nlog2][-t pctcull] [bsdf.sir ..] > bsdf.xml\n",
435 >        "Usage: %s [-g Nlog2][-t pctcull][-l maxlobes] [bsdf.sir ..] > bsdf.xml\n",
436                                  progname);
437          fprintf(stderr,
438          "   or: %s -t{3|4} [-g Nlog2][-t pctcull][{+|-}for[ward]][{+|-}b[ackward]][-e expr][-f file] bsdf_func > bsdf.xml\n",

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines