ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/util/rfluxmtx.c
Revision: 2.45
Committed: Sat May 13 01:50:54 2017 UTC (6 years, 10 months ago) by greg
Content type: text/plain
Branch: MAIN
CVS Tags: rad5R1
Changes since 2.44: +4 -2 lines
Log Message:
Fixed axis reversal in Reinhart sample generator (wrong since August 2014)

File Contents

# User Rev Content
1 greg 2.1 #ifndef lint
2 greg 2.45 static const char RCSid[] = "$Id: rfluxmtx.c,v 2.44 2017/04/13 22:47:17 greg Exp $";
3 greg 2.1 #endif
4     /*
5     * Calculate flux transfer matrix or matrices using rcontrib
6     */
7    
8     #include "copyright.h"
9    
10     #include <ctype.h>
11     #include <stdlib.h>
12     #include "rtio.h"
13     #include "rtmath.h"
14 greg 2.34 #include "paths.h"
15 greg 2.1 #include "bsdf.h"
16     #include "bsdf_m.h"
17     #include "random.h"
18     #include "triangulate.h"
19     #include "platform.h"
20    
21     #define MAXRCARG 512
22    
23     char *progname; /* global argv[0] */
24    
25 greg 2.16 int verbose = 0; /* verbose mode (< 0 no warnings) */
26 greg 2.1
27     char *rcarg[MAXRCARG+1] = {"rcontrib", "-fo+"};
28     int nrcargs = 2;
29    
30     const char overflowerr[] = "%s: too many arguments!\n";
31    
32     #define CHECKARGC(n) if (nrcargs >= MAXRCARG-(n)) \
33     { fprintf(stderr, overflowerr, progname); exit(1); }
34    
35 greg 2.2 int sampcnt = 0; /* sample count (0==unset) */
36 greg 2.1
37 greg 2.2 char *reinhfn = "reinhartb.cal";
38     char *shirchiufn = "disk2square.cal";
39 greg 2.1 char *kfullfn = "klems_full.cal";
40     char *khalffn = "klems_half.cal";
41     char *kquarterfn = "klems_quarter.cal";
42    
43 greg 2.5 /* string indicating parameters */
44     const char PARAMSTART[] = "@rfluxmtx";
45 greg 2.1
46     /* surface type IDs */
47     #define ST_NONE 0
48     #define ST_POLY 1
49     #define ST_RING 2
50     #define ST_SOURCE 3
51    
52     typedef struct surf_s {
53     struct surf_s *next; /* next surface in list */
54     void *priv; /* private data (malloc'ed) */
55     char sname[32]; /* surface name */
56     FVECT snrm; /* surface normal */
57 greg 2.4 double area; /* surface area / proj. solid angle */
58 greg 2.1 short styp; /* surface type */
59 greg 2.2 short nfargs; /* number of real arguments */
60     double farg[1]; /* real values (extends struct) */
61 greg 2.1 } SURF; /* surface structure */
62    
63     typedef struct {
64     FVECT uva[2]; /* tangent axes */
65     int ntris; /* number of triangles */
66     struct ptri {
67 greg 2.41 double afrac; /* fraction of total area */
68 greg 2.1 short vndx[3]; /* vertex indices */
69     } tri[1]; /* triangle array (extends struct) */
70     } POLYTRIS; /* triangulated polygon */
71    
72     typedef struct param_s {
73 greg 2.26 char sign; /* '-' for axis reversal */
74     char hemis[31]; /* hemispherical sampling spec. */
75 greg 2.1 int hsiz; /* hemisphere basis size */
76     int nsurfs; /* number of surfaces */
77     SURF *slist; /* list of surfaces */
78     FVECT vup; /* up vector (zero if unset) */
79     FVECT nrm; /* average normal direction */
80 greg 2.26 FVECT udir, vdir; /* tangent axes */
81 greg 2.1 char *outfn; /* output file name (receiver) */
82     int (*sample_basis)(struct param_s *p, int, FILE *);
83     } PARAMS; /* sender/receiver parameters */
84    
85     PARAMS curparams;
86     char curmod[128];
87 greg 2.6 char newparams[1024];
88 greg 2.1
89     typedef int SURFSAMP(FVECT, SURF *, double);
90    
91     static SURFSAMP ssamp_bad, ssamp_poly, ssamp_ring;
92    
93     SURFSAMP *orig_in_surf[4] = {
94     ssamp_bad, ssamp_poly, ssamp_ring, ssamp_bad
95     };
96    
97     /* Clear parameter set */
98     static void
99     clear_params(PARAMS *p, int reset_only)
100     {
101     while (p->slist != NULL) {
102     SURF *sdel = p->slist;
103     p->slist = sdel->next;
104     if (sdel->priv != NULL)
105     free(sdel->priv);
106     free(sdel);
107     }
108     if (reset_only) {
109     p->nsurfs = 0;
110     memset(p->nrm, 0, sizeof(FVECT));
111     memset(p->vup, 0, sizeof(FVECT));
112     p->outfn = NULL;
113     return;
114     }
115 greg 2.6 memset(p, 0, sizeof(PARAMS));
116 greg 2.1 }
117    
118     /* Get surface type from name */
119     static int
120     surf_type(const char *otype)
121     {
122     if (!strcmp(otype, "polygon"))
123     return(ST_POLY);
124     if (!strcmp(otype, "ring"))
125     return(ST_RING);
126     if (!strcmp(otype, "source"))
127     return(ST_SOURCE);
128     return(ST_NONE);
129     }
130    
131     /* Add arguments to oconv command */
132     static char *
133     oconv_command(int ac, char *av[])
134     {
135 greg 2.13 static char oconvbuf[2048] = "!oconv -f ";
136     char *cp = oconvbuf + 10;
137     char *recv = *av++;
138    
139     if (ac-- <= 0)
140     return(NULL);
141 greg 2.43 if (verbose < 0) { /* turn off warnings */
142     strcpy(cp, "-w- ");
143     cp += 4;
144     }
145 greg 2.1 while (ac-- > 0) {
146 greg 2.13 strcpy(cp, *av++);
147     while (*cp) cp++;
148     *cp++ = ' ';
149 greg 2.30 if (cp >= oconvbuf+(sizeof(oconvbuf)-32))
150     goto overrun;
151     }
152     /* receiver goes last */
153     if (matchany(recv, SPECIALS)) {
154     *cp++ = QUOTCHAR;
155     while (*recv) {
156     if (cp >= oconvbuf+(sizeof(oconvbuf)-3))
157     goto overrun;
158     *cp++ = *recv++;
159 greg 2.1 }
160 greg 2.30 *cp++ = QUOTCHAR;
161     *cp = '\0';
162     } else
163     strcpy(cp, recv);
164 greg 2.1 return(oconvbuf);
165 greg 2.30 overrun:
166     fputs(progname, stderr);
167     fputs(": too many file arguments!\n", stderr);
168     exit(1);
169 greg 2.1 }
170    
171     /* Open a pipe to/from a command given as an argument list */
172     static FILE *
173     popen_arglist(char *av[], char *mode)
174     {
175     char cmd[10240];
176    
177     if (!convert_commandline(cmd, sizeof(cmd), av)) {
178     fputs(progname, stderr);
179     fputs(": command line too long in popen_arglist()\n", stderr);
180     return(NULL);
181     }
182 greg 2.16 if (verbose > 0)
183 greg 2.1 fprintf(stderr, "%s: opening pipe %s: %s\n",
184     progname, (*mode=='w') ? "to" : "from", cmd);
185     return(popen(cmd, mode));
186     }
187    
188 schorsch 2.36 #if defined(_WIN32) || defined(_WIN64)
189 greg 2.1 /* Execute system command (Windows version) */
190     static int
191     my_exec(char *av[])
192     {
193     char cmd[10240];
194    
195     if (!convert_commandline(cmd, sizeof(cmd), av)) {
196     fputs(progname, stderr);
197     fputs(": command line too long in my_exec()\n", stderr);
198     return(1);
199     }
200 greg 2.16 if (verbose > 0)
201 greg 2.1 fprintf(stderr, "%s: running: %s\n", progname, cmd);
202     return(system(cmd));
203     }
204     #else
205     /* Execute system command in our stead (Unix version) */
206     static int
207     my_exec(char *av[])
208     {
209     char *compath;
210    
211     if ((compath = getpath((char *)av[0], getenv("PATH"), X_OK)) == NULL) {
212     fprintf(stderr, "%s: cannot locate %s\n", progname, av[0]);
213     return(1);
214     }
215 greg 2.16 if (verbose > 0) {
216 greg 2.1 char cmd[4096];
217     if (!convert_commandline(cmd, sizeof(cmd), av))
218     strcpy(cmd, "COMMAND TOO LONG TO SHOW");
219     fprintf(stderr, "%s: running: %s\n", progname, cmd);
220     }
221     execv(compath, av); /* successful call never returns */
222     perror(compath);
223     return(1);
224     }
225     #endif
226    
227     /* Get normalized direction vector from string specification */
228     static int
229     get_direction(FVECT dv, const char *s)
230     {
231     int sign = 1;
232     int axis = 0;
233    
234     memset(dv, 0, sizeof(FVECT));
235     nextchar:
236     switch (*s) {
237     case '+':
238     ++s;
239     goto nextchar;
240     case '-':
241     sign = -sign;
242     ++s;
243     goto nextchar;
244     case 'z':
245     case 'Z':
246     ++axis;
247     case 'y':
248     case 'Y':
249     ++axis;
250     case 'x':
251     case 'X':
252     dv[axis] = sign;
253     return(!s[1] | isspace(s[1]));
254     default:
255     break;
256     }
257     #ifdef SMLFLT
258     if (sscanf(s, "%f,%f,%f", &dv[0], &dv[1], &dv[2]) != 3)
259     #else
260     if (sscanf(s, "%lf,%lf,%lf", &dv[0], &dv[1], &dv[2]) != 3)
261     #endif
262     return(0);
263     dv[0] *= (RREAL)sign;
264     return(normalize(dv) > 0);
265     }
266    
267     /* Parse program parameters (directives) */
268     static int
269 greg 2.6 parse_params(PARAMS *p, char *pargs)
270 greg 2.1 {
271     char *cp = pargs;
272     int nparams = 0;
273 greg 2.32 int quot;
274 greg 2.1 int i;
275    
276 greg 2.15 for ( ; ; ) {
277 greg 2.1 switch (*cp++) {
278     case 'h':
279     if (*cp++ != '=')
280     break;
281 greg 2.26 if ((*cp == '+') | (*cp == '-'))
282     p->sign = *cp++;
283     else
284     p->sign = '+';
285 greg 2.6 p->hsiz = 0;
286 greg 2.1 i = 0;
287     while (*cp && !isspace(*cp)) {
288     if (isdigit(*cp))
289 greg 2.6 p->hsiz = 10*p->hsiz + *cp - '0';
290     p->hemis[i++] = *cp++;
291 greg 2.1 }
292     if (!i)
293     break;
294 greg 2.6 p->hemis[i] = '\0';
295     p->hsiz += !p->hsiz;
296 greg 2.1 ++nparams;
297     continue;
298     case 'u':
299     if (*cp++ != '=')
300     break;
301 greg 2.6 if (!get_direction(p->vup, cp))
302 greg 2.1 break;
303 greg 2.15 while (*cp && !isspace(*cp++))
304     ;
305 greg 2.1 ++nparams;
306     continue;
307     case 'o':
308     if (*cp++ != '=')
309     break;
310 greg 2.32 quot = 0;
311     if ((*cp == '"') | (*cp == '\''))
312     quot = *cp++;
313 greg 2.1 i = 0;
314 greg 2.32 while (*cp && (quot ? (*cp != quot) : !isspace(*cp))) {
315     i++; cp++;
316     }
317 greg 2.1 if (!i)
318     break;
319 greg 2.32 if (!*cp) {
320     if (quot)
321     break;
322     cp[1] = '\0';
323     }
324     *cp = '\0';
325 greg 2.6 p->outfn = savqstr(cp-i);
326 greg 2.32 *cp++ = quot ? quot : ' ';
327 greg 2.1 ++nparams;
328     continue;
329     case ' ':
330     case '\t':
331     case '\r':
332 greg 2.21 case '\n':
333 greg 2.15 continue;
334 greg 2.1 case '\0':
335     return(nparams);
336     default:
337     break;
338     }
339 greg 2.15 break;
340     }
341     fprintf(stderr, "%s: bad parameter string: %s", progname, pargs);
342 greg 2.1 exit(1);
343     return(-1); /* pro forma return */
344     }
345    
346     /* Add receiver arguments (directives) corresponding to the current modifier */
347     static void
348     finish_receiver(void)
349     {
350 greg 2.4 char sbuf[256];
351     int uniform = 0;
352 greg 2.2 char *calfn = NULL;
353 greg 2.3 char *params = NULL;
354 greg 2.2 char *binv = NULL;
355     char *binf = NULL;
356     char *nbins = NULL;
357 greg 2.1
358     if (!curmod[0]) {
359     fputs(progname, stderr);
360     fputs(": missing receiver surface!\n", stderr);
361     exit(1);
362     }
363     if (curparams.outfn != NULL) { /* add output file spec. */
364     CHECKARGC(2);
365     rcarg[nrcargs++] = "-o";
366     rcarg[nrcargs++] = curparams.outfn;
367     }
368 greg 2.4 /* check arguments */
369 greg 2.1 if (!curparams.hemis[0]) {
370     fputs(progname, stderr);
371     fputs(": missing hemisphere sampling type!\n", stderr);
372     exit(1);
373     }
374     if (normalize(curparams.nrm) == 0) {
375     fputs(progname, stderr);
376     fputs(": undefined normal for hemisphere sampling\n", stderr);
377     exit(1);
378     }
379 greg 2.10 if (normalize(curparams.vup) == 0) {
380 greg 2.1 if (fabs(curparams.nrm[2]) < .7)
381     curparams.vup[2] = 1;
382     else
383     curparams.vup[1] = 1;
384 greg 2.10 }
385 greg 2.4 /* determine sample type/bin */
386 greg 2.44 if ((tolower(curparams.hemis[0]) == 'u') | (curparams.hemis[0] == '1')) {
387 greg 2.20 sprintf(sbuf, "if(-Dx*%g-Dy*%g-Dz*%g,0,-1)",
388     curparams.nrm[0], curparams.nrm[1], curparams.nrm[2]);
389     binv = savqstr(sbuf);
390     nbins = "1"; /* uniform sampling -- one bin */
391 greg 2.4 uniform = 1;
392 greg 2.1 } else if (tolower(curparams.hemis[0]) == 's' &&
393     tolower(curparams.hemis[1]) == 'c') {
394 greg 2.2 /* assign parameters */
395 greg 2.1 if (curparams.hsiz <= 1) {
396     fputs(progname, stderr);
397     fputs(": missing size for Shirley-Chiu sampling!\n", stderr);
398     exit(1);
399     }
400 greg 2.3 calfn = shirchiufn; shirchiufn = NULL;
401 greg 2.26 sprintf(sbuf, "SCdim=%d,rNx=%g,rNy=%g,rNz=%g,Ux=%g,Uy=%g,Uz=%g,RHS=%c1",
402 greg 2.2 curparams.hsiz,
403     curparams.nrm[0], curparams.nrm[1], curparams.nrm[2],
404 greg 2.26 curparams.vup[0], curparams.vup[1], curparams.vup[2],
405     curparams.sign);
406 greg 2.3 params = savqstr(sbuf);
407 greg 2.2 binv = "scbin";
408     nbins = "SCdim*SCdim";
409 greg 2.1 } else if ((tolower(curparams.hemis[0]) == 'r') |
410     (tolower(curparams.hemis[0]) == 't')) {
411 greg 2.2 calfn = reinhfn; reinhfn = NULL;
412 greg 2.26 sprintf(sbuf, "MF=%d,rNx=%g,rNy=%g,rNz=%g,Ux=%g,Uy=%g,Uz=%g,RHS=%c1",
413 greg 2.3 curparams.hsiz,
414     curparams.nrm[0], curparams.nrm[1], curparams.nrm[2],
415 greg 2.26 curparams.vup[0], curparams.vup[1], curparams.vup[2],
416     curparams.sign);
417 greg 2.3 params = savqstr(sbuf);
418     binv = "rbin";
419     nbins = "Nrbins";
420 greg 2.1 } else if (tolower(curparams.hemis[0]) == 'k' &&
421     !curparams.hemis[1] |
422     (tolower(curparams.hemis[1]) == 'f') |
423     (curparams.hemis[1] == '1')) {
424 greg 2.2 calfn = kfullfn; kfullfn = NULL;
425     binf = "kbin";
426     nbins = "Nkbins";
427 greg 2.1 } else if (tolower(curparams.hemis[0]) == 'k' &&
428     (tolower(curparams.hemis[1]) == 'h') |
429     (curparams.hemis[1] == '2')) {
430 greg 2.2 calfn = khalffn; khalffn = NULL;
431     binf = "khbin";
432     nbins = "Nkhbins";
433 greg 2.1 } else if (tolower(curparams.hemis[0]) == 'k' &&
434     (tolower(curparams.hemis[1]) == 'q') |
435     (curparams.hemis[1] == '4')) {
436 greg 2.2 calfn = kquarterfn; kquarterfn = NULL;
437     binf = "kqbin";
438     nbins = "Nkqbins";
439 greg 2.1 } else {
440     fprintf(stderr, "%s: unrecognized hemisphere sampling: h=%s\n",
441     progname, curparams.hemis);
442     exit(1);
443     }
444 greg 2.26 if (tolower(curparams.hemis[0]) == 'k') {
445     sprintf(sbuf, "RHS=%c1", curparams.sign);
446     params = savqstr(sbuf);
447     }
448 greg 2.4 if (!uniform & (curparams.slist->styp == ST_SOURCE)) {
449     SURF *sp;
450     for (sp = curparams.slist; sp != NULL; sp = sp->next)
451     if (fabs(sp->area - PI) > 1e-3) {
452     fprintf(stderr, "%s: source '%s' must be 180-degrees\n",
453     progname, sp->sname);
454     exit(1);
455     }
456     }
457 greg 2.2 if (calfn != NULL) { /* add cal file if needed */
458     CHECKARGC(2);
459     rcarg[nrcargs++] = "-f";
460     rcarg[nrcargs++] = calfn;
461     }
462 greg 2.3 if (params != NULL) { /* parameters _after_ cal file */
463     CHECKARGC(2);
464     rcarg[nrcargs++] = "-p";
465     rcarg[nrcargs++] = params;
466     }
467 greg 2.2 if (nbins != NULL) { /* add #bins if set */
468     CHECKARGC(2);
469     rcarg[nrcargs++] = "-bn";
470     rcarg[nrcargs++] = nbins;
471     }
472 greg 2.3 if (binv != NULL) {
473 greg 2.2 CHECKARGC(2); /* assign bin variable */
474     rcarg[nrcargs++] = "-b";
475     rcarg[nrcargs++] = binv;
476     } else if (binf != NULL) {
477     CHECKARGC(2); /* assign bin function */
478 greg 2.3 rcarg[nrcargs++] = "-b";
479 greg 2.2 sprintf(sbuf, "%s(%g,%g,%g,%g,%g,%g)", binf,
480     curparams.nrm[0], curparams.nrm[1], curparams.nrm[2],
481     curparams.vup[0], curparams.vup[1], curparams.vup[2]);
482     rcarg[nrcargs++] = savqstr(sbuf);
483     }
484     CHECKARGC(2); /* modifier argument goes last */
485 greg 2.1 rcarg[nrcargs++] = "-m";
486     rcarg[nrcargs++] = savqstr(curmod);
487     }
488    
489     /* Make randomly oriented tangent plane axes for given normal direction */
490     static void
491     make_axes(FVECT uva[2], const FVECT nrm)
492     {
493     int i;
494    
495 greg 2.29 if (!getperpendicular(uva[0], nrm, 1)) {
496 greg 2.1 fputs(progname, stderr);
497     fputs(": bad surface normal in make_axes!\n", stderr);
498     exit(1);
499     }
500 greg 2.19 fcross(uva[1], nrm, uva[0]);
501 greg 2.1 }
502    
503     /* Illegal sender surfaces end up here */
504     static int
505     ssamp_bad(FVECT orig, SURF *sp, double x)
506     {
507     fprintf(stderr, "%s: illegal sender surface '%s'\n",
508     progname, sp->sname);
509     return(0);
510     }
511    
512     /* Generate origin on ring surface from uniform random variable */
513     static int
514     ssamp_ring(FVECT orig, SURF *sp, double x)
515     {
516     FVECT *uva = (FVECT *)sp->priv;
517     double samp2[2];
518     double uv[2];
519     int i;
520    
521     if (uva == NULL) { /* need tangent axes */
522     uva = (FVECT *)malloc(sizeof(FVECT)*2);
523     if (uva == NULL) {
524     fputs(progname, stderr);
525     fputs(": out of memory in ssamp_ring!\n", stderr);
526     return(0);
527     }
528     make_axes(uva, sp->snrm);
529     sp->priv = (void *)uva;
530     }
531     SDmultiSamp(samp2, 2, x);
532 greg 2.6 samp2[0] = sqrt(samp2[0]*sp->area*(1./PI) + sp->farg[6]*sp->farg[6]);
533 greg 2.1 samp2[1] *= 2.*PI;
534     uv[0] = samp2[0]*tcos(samp2[1]);
535     uv[1] = samp2[0]*tsin(samp2[1]);
536     for (i = 3; i--; )
537     orig[i] = sp->farg[i] + uv[0]*uva[0][i] + uv[1]*uva[1][i];
538     return(1);
539     }
540    
541     /* Add triangle to polygon's list (call-back function) */
542     static int
543     add_triangle(const Vert2_list *tp, int a, int b, int c)
544     {
545     POLYTRIS *ptp = (POLYTRIS *)tp->p;
546     struct ptri *trip = ptp->tri + ptp->ntris++;
547    
548     trip->vndx[0] = a;
549     trip->vndx[1] = b;
550     trip->vndx[2] = c;
551     return(1);
552     }
553    
554     /* Generate origin on polygon surface from uniform random variable */
555     static int
556     ssamp_poly(FVECT orig, SURF *sp, double x)
557     {
558     POLYTRIS *ptp = (POLYTRIS *)sp->priv;
559     double samp2[2];
560     double *v0, *v1, *v2;
561     int i;
562    
563     if (ptp == NULL) { /* need to triangulate */
564     ptp = (POLYTRIS *)malloc(sizeof(POLYTRIS) +
565 greg 2.40 sizeof(struct ptri)*(sp->nfargs/3 - 3));
566 greg 2.1 if (ptp == NULL)
567     goto memerr;
568     if (sp->nfargs == 3) { /* simple case */
569     ptp->ntris = 1;
570     ptp->tri[0].vndx[0] = 0;
571     ptp->tri[0].vndx[1] = 1;
572     ptp->tri[0].vndx[2] = 2;
573     ptp->tri[0].afrac = 1;
574     } else {
575     Vert2_list *v2l = polyAlloc(sp->nfargs/3);
576     if (v2l == NULL)
577     goto memerr;
578     make_axes(ptp->uva, sp->snrm);
579     for (i = v2l->nv; i--; ) {
580     v2l->v[i].mX = DOT(sp->farg+3*i, ptp->uva[0]);
581     v2l->v[i].mY = DOT(sp->farg+3*i, ptp->uva[1]);
582     }
583     ptp->ntris = 0;
584     v2l->p = (void *)ptp;
585 greg 2.7 if (!polyTriangulate(v2l, add_triangle)) {
586     fprintf(stderr,
587     "%s: cannot triangulate polygon '%s'\n",
588     progname, sp->sname);
589 greg 2.1 return(0);
590 greg 2.7 }
591 greg 2.1 for (i = ptp->ntris; i--; ) {
592     int a = ptp->tri[i].vndx[0];
593     int b = ptp->tri[i].vndx[1];
594     int c = ptp->tri[i].vndx[2];
595     ptp->tri[i].afrac =
596     (v2l->v[a].mX*v2l->v[b].mY -
597     v2l->v[b].mX*v2l->v[a].mY +
598     v2l->v[b].mX*v2l->v[c].mY -
599     v2l->v[c].mX*v2l->v[b].mY +
600     v2l->v[c].mX*v2l->v[a].mY -
601     v2l->v[a].mX*v2l->v[c].mY) /
602     (2.*sp->area);
603     }
604     polyFree(v2l);
605     }
606     sp->priv = (void *)ptp;
607     }
608     /* pick triangle by partial area */
609 greg 2.40 for (i = 0; i < ptp->ntris-1 && x > ptp->tri[i].afrac; i++)
610 greg 2.1 x -= ptp->tri[i].afrac;
611     SDmultiSamp(samp2, 2, x/ptp->tri[i].afrac);
612     samp2[0] *= samp2[1] = sqrt(samp2[1]);
613     samp2[1] = 1. - samp2[1];
614     v0 = sp->farg + 3*ptp->tri[i].vndx[0];
615     v1 = sp->farg + 3*ptp->tri[i].vndx[1];
616     v2 = sp->farg + 3*ptp->tri[i].vndx[2];
617     for (i = 3; i--; )
618     orig[i] = v0[i] + samp2[0]*(v1[i] - v0[i])
619     + samp2[1]*(v2[i] - v0[i]) ;
620     return(1);
621     memerr:
622     fputs(progname, stderr);
623     fputs(": out of memory in ssamp_poly!\n", stderr);
624     return(0);
625     }
626    
627     /* Compute sample origin based on projected areas of sender subsurfaces */
628     static int
629     sample_origin(PARAMS *p, FVECT orig, const FVECT rdir, double x)
630     {
631     static double *projsa;
632     static int nall;
633     double tarea = 0;
634     int i;
635     SURF *sp;
636     /* special case for lone surface */
637     if (p->nsurfs == 1) {
638     sp = p->slist;
639 greg 2.31 if (DOT(sp->snrm, rdir) >= FTINY) {
640 greg 2.7 fprintf(stderr,
641     "%s: internal - sample behind sender '%s'\n",
642     progname, sp->sname);
643     return(0);
644     }
645 greg 2.1 return((*orig_in_surf[sp->styp])(orig, sp, x));
646     }
647     if (p->nsurfs > nall) { /* (re)allocate surface area cache */
648     if (projsa) free(projsa);
649     projsa = (double *)malloc(sizeof(double)*p->nsurfs);
650 greg 2.37 if (projsa == NULL) {
651     fputs(progname, stderr);
652     fputs(": out of memory in sample_origin!\n", stderr);
653     exit(1);
654     }
655 greg 2.1 nall = p->nsurfs;
656     }
657     /* compute projected areas */
658     for (i = 0, sp = p->slist; sp != NULL; i++, sp = sp->next) {
659     projsa[i] = -DOT(sp->snrm, rdir) * sp->area;
660     tarea += projsa[i] *= (double)(projsa[i] > FTINY);
661     }
662 greg 2.7 if (tarea <= FTINY) { /* wrong side of sender? */
663     fputs(progname, stderr);
664     fputs(": internal - sample behind all sender elements!\n",
665     stderr);
666 greg 2.1 return(0);
667 greg 2.7 }
668 greg 2.1 tarea *= x; /* get surface from list */
669     for (i = 0, sp = p->slist; tarea > projsa[i]; sp = sp->next)
670     tarea -= projsa[i++];
671     return((*orig_in_surf[sp->styp])(orig, sp, tarea/projsa[i]));
672     }
673    
674     /* Uniform sample generator */
675     static int
676     sample_uniform(PARAMS *p, int b, FILE *fp)
677     {
678     int n = sampcnt;
679     double samp3[3];
680     double duvw[3];
681     FVECT orig_dir[2];
682     int i;
683    
684     if (fp == NULL) /* just requesting number of bins? */
685     return(1);
686    
687     while (n--) { /* stratified hemisphere sampling */
688     SDmultiSamp(samp3, 3, (n+frandom())/sampcnt);
689     SDsquare2disk(duvw, samp3[1], samp3[2]);
690     duvw[2] = -sqrt(1. - duvw[0]*duvw[0] - duvw[1]*duvw[1]);
691     for (i = 3; i--; )
692     orig_dir[1][i] = duvw[0]*p->udir[i] +
693     duvw[1]*p->vdir[i] +
694     duvw[2]*p->nrm[i] ;
695     if (!sample_origin(p, orig_dir[0], orig_dir[1], samp3[0]))
696     return(0);
697 greg 2.38 if (putbinary(orig_dir, sizeof(FVECT), 2, fp) != 2)
698 greg 2.1 return(0);
699     }
700     return(1);
701     }
702    
703     /* Shirly-Chiu sample generator */
704     static int
705     sample_shirchiu(PARAMS *p, int b, FILE *fp)
706     {
707     int n = sampcnt;
708     double samp3[3];
709     double duvw[3];
710     FVECT orig_dir[2];
711     int i;
712    
713     if (fp == NULL) /* just requesting number of bins? */
714     return(p->hsiz*p->hsiz);
715    
716     while (n--) { /* stratified sampling */
717     SDmultiSamp(samp3, 3, (n+frandom())/sampcnt);
718     SDsquare2disk(duvw, (b/p->hsiz + samp3[1])/curparams.hsiz,
719     (b%p->hsiz + samp3[2])/curparams.hsiz);
720     duvw[2] = sqrt(1. - duvw[0]*duvw[0] - duvw[1]*duvw[1]);
721     for (i = 3; i--; )
722     orig_dir[1][i] = -duvw[0]*p->udir[i] -
723     duvw[1]*p->vdir[i] -
724     duvw[2]*p->nrm[i] ;
725     if (!sample_origin(p, orig_dir[0], orig_dir[1], samp3[0]))
726     return(0);
727 greg 2.38 if (putbinary(orig_dir, sizeof(FVECT), 2, fp) != 2)
728 greg 2.1 return(0);
729     }
730     return(1);
731     }
732    
733     /* Reinhart/Tregenza sample generator */
734     static int
735     sample_reinhart(PARAMS *p, int b, FILE *fp)
736     {
737     #define T_NALT 7
738     static const int tnaz[T_NALT] = {30, 30, 24, 24, 18, 12, 6};
739     const int RowMax = T_NALT*p->hsiz + 1;
740 greg 2.7 const double RAH = (.5*PI)/(RowMax-.5);
741 greg 2.1 #define rnaz(r) (r >= RowMax-1 ? 1 : p->hsiz*tnaz[r/p->hsiz])
742     int n = sampcnt;
743     int row, col;
744     double samp3[3];
745     double alt, azi;
746     double duvw[3];
747     FVECT orig_dir[2];
748     int i;
749    
750     if (fp == NULL) { /* just requesting number of bins? */
751     n = 0;
752     for (row = RowMax; row--; ) n += rnaz(row);
753     return(n);
754     }
755     row = 0; /* identify row & column */
756     col = b;
757     while (col >= rnaz(row)) {
758     col -= rnaz(row);
759     ++row;
760     }
761     while (n--) { /* stratified sampling */
762     SDmultiSamp(samp3, 3, (n+frandom())/sampcnt);
763 greg 2.45 if (row >= RowMax-1) /* avoid crowding at zenith */
764     samp3[1] *= samp3[1];
765 greg 2.1 alt = (row+samp3[1])*RAH;
766     azi = (2.*PI)*(col+samp3[2]-.5)/rnaz(row);
767 greg 2.7 duvw[2] = cos(alt); /* measured from horizon */
768 greg 2.18 duvw[0] = tsin(azi)*duvw[2];
769 greg 2.45 duvw[1] = -tcos(azi)*duvw[2];
770 greg 2.1 duvw[2] = sqrt(1. - duvw[2]*duvw[2]);
771     for (i = 3; i--; )
772     orig_dir[1][i] = -duvw[0]*p->udir[i] -
773     duvw[1]*p->vdir[i] -
774     duvw[2]*p->nrm[i] ;
775     if (!sample_origin(p, orig_dir[0], orig_dir[1], samp3[0]))
776     return(0);
777 greg 2.38 if (putbinary(orig_dir, sizeof(FVECT), 2, fp) != 2)
778 greg 2.1 return(0);
779     }
780     return(1);
781     #undef rnaz
782     #undef T_NALT
783     }
784    
785     /* Klems sample generator */
786     static int
787     sample_klems(PARAMS *p, int b, FILE *fp)
788     {
789     static const char bname[4][20] = {
790     "LBNL/Klems Full",
791     "LBNL/Klems Half",
792     "INTERNAL ERROR",
793     "LBNL/Klems Quarter"
794     };
795     static ANGLE_BASIS *kbasis[4];
796     const int bi = p->hemis[1] - '1';
797     int n = sampcnt;
798     double samp2[2];
799     double duvw[3];
800     FVECT orig_dir[2];
801     int i;
802    
803     if (!kbasis[bi]) { /* need to get basis, first */
804     for (i = 4; i--; )
805     if (!strcasecmp(abase_list[i].name, bname[bi])) {
806     kbasis[bi] = &abase_list[i];
807     break;
808     }
809     if (i < 0) {
810     fprintf(stderr, "%s: unknown hemisphere basis '%s'\n",
811     progname, bname[bi]);
812     return(0);
813     }
814     }
815     if (fp == NULL) /* just requesting number of bins? */
816     return(kbasis[bi]->nangles);
817    
818     while (n--) { /* stratified sampling */
819     SDmultiSamp(samp2, 2, (n+frandom())/sampcnt);
820 greg 2.27 if (!fo_getvec(duvw, b+samp2[1], kbasis[bi]))
821 greg 2.1 return(0);
822     for (i = 3; i--; )
823 greg 2.26 orig_dir[1][i] = -duvw[0]*p->udir[i] -
824     duvw[1]*p->vdir[i] -
825 greg 2.1 duvw[2]*p->nrm[i] ;
826     if (!sample_origin(p, orig_dir[0], orig_dir[1], samp2[0]))
827     return(0);
828 greg 2.38 if (putbinary(orig_dir, sizeof(FVECT), 2, fp) != 2)
829 greg 2.1 return(0);
830     }
831     return(1);
832     }
833    
834     /* Prepare hemisphere basis sampler that will send rays to rcontrib */
835     static int
836     prepare_sampler(void)
837     {
838     if (curparams.slist == NULL) { /* missing sample surface! */
839     fputs(progname, stderr);
840     fputs(": no sender surface!\n", stderr);
841     return(-1);
842     }
843 greg 2.16 /* misplaced output file spec. */
844     if ((curparams.outfn != NULL) & (verbose >= 0))
845 greg 2.1 fprintf(stderr, "%s: warning - ignoring output file in sender ('%s')\n",
846     progname, curparams.outfn);
847     /* check/set basis hemisphere */
848     if (!curparams.hemis[0]) {
849     fputs(progname, stderr);
850     fputs(": missing sender sampling type!\n", stderr);
851     return(-1);
852     }
853     if (normalize(curparams.nrm) == 0) {
854     fputs(progname, stderr);
855     fputs(": undefined normal for sender sampling\n", stderr);
856     return(-1);
857     }
858 greg 2.10 if (normalize(curparams.vup) == 0) {
859 greg 2.1 if (fabs(curparams.nrm[2]) < .7)
860     curparams.vup[2] = 1;
861     else
862     curparams.vup[1] = 1;
863 greg 2.10 }
864 greg 2.26 fcross(curparams.udir, curparams.vup, curparams.nrm);
865 greg 2.1 if (normalize(curparams.udir) == 0) {
866     fputs(progname, stderr);
867     fputs(": up vector coincides with sender normal\n", stderr);
868     return(-1);
869     }
870 greg 2.26 fcross(curparams.vdir, curparams.nrm, curparams.udir);
871     if (curparams.sign == '-') { /* left-handed coordinate system? */
872     curparams.udir[0] *= -1.;
873     curparams.udir[1] *= -1.;
874     curparams.udir[2] *= -1.;
875     }
876 greg 2.44 if ((tolower(curparams.hemis[0]) == 'u') | (curparams.hemis[0] == '1'))
877 greg 2.1 curparams.sample_basis = sample_uniform;
878     else if (tolower(curparams.hemis[0]) == 's' &&
879     tolower(curparams.hemis[1]) == 'c')
880     curparams.sample_basis = sample_shirchiu;
881     else if ((tolower(curparams.hemis[0]) == 'r') |
882     (tolower(curparams.hemis[0]) == 't'))
883     curparams.sample_basis = sample_reinhart;
884     else if (tolower(curparams.hemis[0]) == 'k') {
885     switch (curparams.hemis[1]) {
886     case '1':
887     case '2':
888     case '4':
889     break;
890     case 'f':
891     case 'F':
892     case '\0':
893     curparams.hemis[1] = '1';
894     break;
895     case 'h':
896     case 'H':
897     curparams.hemis[1] = '2';
898     break;
899     case 'q':
900     case 'Q':
901     curparams.hemis[1] = '4';
902     break;
903     default:
904     goto unrecognized;
905     }
906     curparams.hemis[2] = '\0';
907     curparams.sample_basis = sample_klems;
908     } else
909     goto unrecognized;
910     /* return number of bins */
911     return((*curparams.sample_basis)(&curparams,0,NULL));
912     unrecognized:
913     fprintf(stderr, "%s: unrecognized sender sampling: h=%s\n",
914     progname, curparams.hemis);
915     return(-1);
916     }
917    
918     /* Compute normal and area for polygon */
919     static int
920     finish_polygon(SURF *p)
921     {
922     const int nv = p->nfargs / 3;
923     FVECT e1, e2, vc;
924     int i;
925    
926     memset(p->snrm, 0, sizeof(FVECT));
927     VSUB(e1, p->farg+3, p->farg);
928     for (i = 2; i < nv; i++) {
929     VSUB(e2, p->farg+3*i, p->farg);
930     VCROSS(vc, e1, e2);
931     p->snrm[0] += vc[0];
932     p->snrm[1] += vc[1];
933     p->snrm[2] += vc[2];
934     VCOPY(e1, e2);
935     }
936     p->area = normalize(p->snrm)*0.5;
937     return(p->area > FTINY);
938     }
939    
940     /* Add a surface to our current parameters */
941     static void
942     add_surface(int st, const char *oname, FILE *fp)
943     {
944     SURF *snew;
945     int n;
946     /* get floating-point arguments */
947     if (!fscanf(fp, "%d", &n)) return;
948     while (n-- > 0) fscanf(fp, "%*s");
949     if (!fscanf(fp, "%d", &n)) return;
950     while (n-- > 0) fscanf(fp, "%*d");
951     if (!fscanf(fp, "%d", &n) || n <= 0) return;
952     snew = (SURF *)malloc(sizeof(SURF) + sizeof(double)*(n-1));
953     if (snew == NULL) {
954     fputs(progname, stderr);
955 greg 2.37 fputs(": out of memory in add_surface!\n", stderr);
956 greg 2.1 exit(1);
957     }
958     strncpy(snew->sname, oname, sizeof(snew->sname)-1);
959     snew->sname[sizeof(snew->sname)-1] = '\0';
960     snew->styp = st;
961     snew->priv = NULL;
962     snew->nfargs = n;
963     for (n = 0; n < snew->nfargs; n++)
964     if (fscanf(fp, "%lf", &snew->farg[n]) != 1) {
965     fprintf(stderr, "%s: error reading arguments for '%s'\n",
966     progname, oname);
967     exit(1);
968     }
969     switch (st) {
970     case ST_RING:
971     if (snew->nfargs != 8)
972     goto badcount;
973     VCOPY(snew->snrm, snew->farg+3);
974     if (normalize(snew->snrm) == 0)
975     goto badnorm;
976     if (snew->farg[7] < snew->farg[6]) {
977     double t = snew->farg[7];
978     snew->farg[7] = snew->farg[6];
979     snew->farg[6] = t;
980     }
981     snew->area = PI*(snew->farg[7]*snew->farg[7] -
982     snew->farg[6]*snew->farg[6]);
983     break;
984     case ST_POLY:
985     if (snew->nfargs < 9 || snew->nfargs % 3)
986     goto badcount;
987     finish_polygon(snew);
988     break;
989     case ST_SOURCE:
990     if (snew->nfargs != 4)
991     goto badcount;
992 greg 2.8 for (n = 3; n--; ) /* need to reverse "normal" */
993     snew->snrm[n] = -snew->farg[n];
994 greg 2.1 if (normalize(snew->snrm) == 0)
995     goto badnorm;
996 greg 2.4 snew->area = sin((PI/180./2.)*snew->farg[3]);
997     snew->area *= PI*snew->area;
998 greg 2.1 break;
999     }
1000 greg 2.16 if ((snew->area <= FTINY) & (verbose >= 0)) {
1001 greg 2.1 fprintf(stderr, "%s: warning - zero area for surface '%s'\n",
1002     progname, oname);
1003     free(snew);
1004     return;
1005     }
1006     VSUM(curparams.nrm, curparams.nrm, snew->snrm, snew->area);
1007     snew->next = curparams.slist;
1008     curparams.slist = snew;
1009     curparams.nsurfs++;
1010     return;
1011     badcount:
1012 greg 2.7 fprintf(stderr, "%s: bad argument count for surface element '%s'\n",
1013 greg 2.1 progname, oname);
1014     exit(1);
1015     badnorm:
1016 greg 2.7 fprintf(stderr, "%s: bad orientation for surface element '%s'\n",
1017 greg 2.1 progname, oname);
1018     exit(1);
1019     }
1020    
1021     /* Parse a receiver object (look for modifiers to add to rcontrib command) */
1022     static int
1023     add_recv_object(FILE *fp)
1024     {
1025     int st;
1026     char thismod[128], otype[32], oname[128];
1027     int n;
1028    
1029     if (fscanf(fp, "%s %s %s", thismod, otype, oname) != 3)
1030     return(0); /* must have hit EOF! */
1031     if (!strcmp(otype, "alias")) {
1032     fscanf(fp, "%*s"); /* skip alias */
1033     return(0);
1034     }
1035     /* is it a new receiver? */
1036     if ((st = surf_type(otype)) != ST_NONE) {
1037     if (curparams.slist != NULL && (st == ST_SOURCE) ^
1038     (curparams.slist->styp == ST_SOURCE)) {
1039     fputs(progname, stderr);
1040     fputs(": cannot mix source/non-source receivers!\n", stderr);
1041     return(-1);
1042     }
1043     if (strcmp(thismod, curmod)) {
1044     if (curmod[0]) { /* output last receiver? */
1045     finish_receiver();
1046     clear_params(&curparams, 1);
1047     }
1048 greg 2.6 parse_params(&curparams, newparams);
1049     newparams[0] = '\0';
1050 greg 2.1 strcpy(curmod, thismod);
1051     }
1052     add_surface(st, oname, fp); /* read & store surface */
1053     return(1);
1054     }
1055     /* else skip arguments */
1056 greg 2.8 if (!fscanf(fp, "%d", &n)) return(0);
1057 greg 2.1 while (n-- > 0) fscanf(fp, "%*s");
1058 greg 2.10 if (!fscanf(fp, "%d", &n)) return(0);
1059 greg 2.1 while (n-- > 0) fscanf(fp, "%*d");
1060 greg 2.10 if (!fscanf(fp, "%d", &n)) return(0);
1061 greg 2.1 while (n-- > 0) fscanf(fp, "%*f");
1062     return(0);
1063     }
1064    
1065     /* Parse a sender object */
1066     static int
1067     add_send_object(FILE *fp)
1068     {
1069     int st;
1070 greg 2.24 char thismod[128], otype[32], oname[128];
1071 greg 2.1 int n;
1072    
1073 greg 2.25 if (fscanf(fp, "%s %s %s", thismod, otype, oname) != 3)
1074 greg 2.1 return(0); /* must have hit EOF! */
1075     if (!strcmp(otype, "alias")) {
1076     fscanf(fp, "%*s"); /* skip alias */
1077     return(0);
1078     }
1079     /* is it a new surface? */
1080     if ((st = surf_type(otype)) != ST_NONE) {
1081     if (st == ST_SOURCE) {
1082     fputs(progname, stderr);
1083     fputs(": cannot use source as a sender!\n", stderr);
1084     return(-1);
1085     }
1086 greg 2.24 if (strcmp(thismod, curmod)) {
1087     if (curmod[0]) {
1088     fputs(progname, stderr);
1089     fputs(": warning - multiple modifiers in sender\n",
1090     stderr);
1091     }
1092     strcpy(curmod, thismod);
1093     }
1094 greg 2.6 parse_params(&curparams, newparams);
1095     newparams[0] = '\0';
1096 greg 2.1 add_surface(st, oname, fp); /* read & store surface */
1097     return(0);
1098     }
1099     /* else skip arguments */
1100 greg 2.8 if (!fscanf(fp, "%d", &n)) return(0);
1101 greg 2.1 while (n-- > 0) fscanf(fp, "%*s");
1102 greg 2.10 if (!fscanf(fp, "%d", &n)) return(0);
1103 greg 2.1 while (n-- > 0) fscanf(fp, "%*d");
1104 greg 2.10 if (!fscanf(fp, "%d", &n)) return(0);
1105 greg 2.1 while (n-- > 0) fscanf(fp, "%*f");
1106     return(0);
1107     }
1108    
1109     /* Load a Radiance scene using the given callback function for objects */
1110     static int
1111     load_scene(const char *inspec, int (*ocb)(FILE *))
1112     {
1113     int rv = 0;
1114     char inpbuf[1024];
1115     FILE *fp;
1116     int c;
1117    
1118     if (*inspec == '!')
1119     fp = popen(inspec+1, "r");
1120     else
1121     fp = fopen(inspec, "r");
1122     if (fp == NULL) {
1123     fprintf(stderr, "%s: cannot load '%s'\n", progname, inspec);
1124     return(-1);
1125     }
1126     while ((c = getc(fp)) != EOF) { /* load receiver data */
1127     if (isspace(c)) /* skip leading white space */
1128     continue;
1129     if (c == '!') { /* read from a new command */
1130     inpbuf[0] = c;
1131     if (fgetline(inpbuf+1, sizeof(inpbuf)-1, fp) != NULL) {
1132     if ((c = load_scene(inpbuf, ocb)) < 0)
1133     return(c);
1134     rv += c;
1135     }
1136     continue;
1137     }
1138     if (c == '#') { /* parameters/comment */
1139 greg 2.5 if ((c = getc(fp)) == EOF || ungetc(c,fp) == EOF)
1140     break;
1141     if (!isspace(c) && fscanf(fp, "%s", inpbuf) == 1 &&
1142 greg 2.1 !strcmp(inpbuf, PARAMSTART)) {
1143     if (fgets(inpbuf, sizeof(inpbuf), fp) != NULL)
1144 greg 2.6 strcat(newparams, inpbuf);
1145 greg 2.1 continue;
1146     }
1147 greg 2.10 while ((c = getc(fp)) != EOF && c != '\n')
1148 greg 2.1 ; /* else skipping comment */
1149     continue;
1150     }
1151     ungetc(c, fp); /* else check object for receiver */
1152     c = (*ocb)(fp);
1153     if (c < 0)
1154     return(c);
1155     rv += c;
1156     }
1157     /* close our input stream */
1158     c = (*inspec == '!') ? pclose(fp) : fclose(fp);
1159     if (c != 0) {
1160     fprintf(stderr, "%s: error loading '%s'\n", progname, inspec);
1161     return(-1);
1162     }
1163     return(rv);
1164     }
1165    
1166     /* Get command arguments and run program */
1167     int
1168     main(int argc, char *argv[])
1169     {
1170     char fmtopt[6] = "-faa"; /* default output is ASCII */
1171 greg 2.4 char *xrs=NULL, *yrs=NULL, *ldopt=NULL;
1172 greg 2.11 char *iropt = NULL;
1173 greg 2.1 char *sendfn;
1174     char sampcntbuf[32], nsbinbuf[32];
1175     FILE *rcfp;
1176     int nsbins;
1177     int a, i;
1178     /* screen rcontrib options */
1179     progname = argv[0];
1180 greg 2.14 for (a = 1; a < argc-2; a++) {
1181     int na;
1182     /* check for argument expansion */
1183     while ((na = expandarg(&argc, &argv, a)) > 0)
1184     ;
1185     if (na < 0) {
1186     fprintf(stderr, "%s: cannot expand '%s'\n",
1187     progname, argv[a]);
1188     return(1);
1189     }
1190     if (argv[a][0] != '-' || !argv[a][1])
1191     break;
1192     na = 1;
1193     switch (argv[a][1]) { /* !! Keep consistent !! */
1194 greg 2.1 case 'v': /* verbose mode */
1195 greg 2.16 verbose = 1;
1196 greg 2.1 na = 0;
1197     continue;
1198     case 'f': /* special case for -fo, -ff, etc. */
1199     switch (argv[a][2]) {
1200     case '\0': /* cal file */
1201     goto userr;
1202     case 'o': /* force output */
1203     goto userr;
1204     case 'a': /* output format */
1205     case 'f':
1206     case 'd':
1207     case 'c':
1208 greg 2.6 if (!(fmtopt[3] = argv[a][3]))
1209     fmtopt[3] = argv[a][2];
1210     fmtopt[2] = argv[a][2];
1211 greg 2.1 na = 0;
1212     continue; /* will pass later */
1213     default:
1214     goto userr;
1215     }
1216     break;
1217 greg 2.4 case 'x': /* x-resolution */
1218     xrs = argv[++a];
1219     na = 0;
1220     continue;
1221     case 'y': /* y-resolution */
1222     yrs = argv[++a];
1223     na = 0;
1224     continue;
1225 greg 2.1 case 'c': /* number of samples */
1226 greg 2.6 sampcnt = atoi(argv[++a]);
1227 greg 2.1 if (sampcnt <= 0)
1228     goto userr;
1229     na = 0; /* we re-add this later */
1230     continue;
1231 greg 2.8 case 'I': /* only for pass-through mode */
1232 greg 2.11 case 'i':
1233     iropt = argv[a];
1234 greg 2.8 na = 0;
1235     continue;
1236 greg 2.16 case 'w': /* options without arguments */
1237     if (argv[a][2] != '+') verbose = -1;
1238     case 'V':
1239 greg 2.1 case 'u':
1240     case 'h':
1241 greg 2.4 case 'r':
1242 greg 2.1 break;
1243     case 'n': /* options with 1 argument */
1244     case 's':
1245     case 'o':
1246     na = 2;
1247     break;
1248     case 'b': /* special case */
1249     if (argv[a][2] != 'v') goto userr;
1250     break;
1251     case 'l': /* special case */
1252 greg 2.4 if (argv[a][2] == 'd') {
1253     ldopt = argv[a];
1254     na = 0;
1255     continue;
1256     }
1257 greg 2.1 na = 2;
1258     break;
1259     case 'd': /* special case */
1260     if (argv[a][2] != 'v') na = 2;
1261     break;
1262     case 'a': /* special case */
1263 greg 2.28 if (argv[a][2] == 'p') {
1264     na = 2; /* photon map [+ bandwidth(s)] */
1265     if (a < argc-3 && atoi(argv[a+1]) > 0)
1266     na += 1 + (a < argc-4 && atoi(argv[a+2]) > 0);
1267     } else
1268     na = (argv[a][2] == 'v') ? 4 : 2;
1269 greg 2.1 break;
1270     case 'm': /* special case */
1271     if (!argv[a][2]) goto userr;
1272     na = (argv[a][2] == 'e') | (argv[a][2] == 'a') ? 4 : 2;
1273     break;
1274     default: /* anything else is verbotten */
1275     goto userr;
1276     }
1277     if (na <= 0) continue;
1278     CHECKARGC(na); /* pass on option */
1279     rcarg[nrcargs++] = argv[a];
1280     while (--na) /* + arguments if any */
1281     rcarg[nrcargs++] = argv[++a];
1282     }
1283     if (a > argc-2)
1284     goto userr; /* check at end of options */
1285     sendfn = argv[a++]; /* assign sender & receiver inputs */
1286     if (sendfn[0] == '-') { /* user wants pass-through mode? */
1287     if (sendfn[1]) goto userr;
1288     sendfn = NULL;
1289 greg 2.11 if (iropt) {
1290 greg 2.8 CHECKARGC(1);
1291 greg 2.11 rcarg[nrcargs++] = iropt;
1292 greg 2.8 }
1293 greg 2.4 if (xrs) {
1294     CHECKARGC(2);
1295     rcarg[nrcargs++] = "-x";
1296     rcarg[nrcargs++] = xrs;
1297     }
1298     if (yrs) {
1299     CHECKARGC(2);
1300     rcarg[nrcargs++] = "-y";
1301     rcarg[nrcargs++] = yrs;
1302     }
1303     if (ldopt) {
1304     CHECKARGC(1);
1305     rcarg[nrcargs++] = ldopt;
1306     }
1307 greg 2.1 if (sampcnt <= 0) sampcnt = 1;
1308 greg 2.8 } else { /* else in sampling mode */
1309 greg 2.11 if (iropt) {
1310 greg 2.8 fputs(progname, stderr);
1311 greg 2.11 fputs(": -i, -I supported for pass-through only\n", stderr);
1312 greg 2.8 return(1);
1313     }
1314 greg 2.42 #ifdef SMLFLT
1315     fmtopt[2] = 'f';
1316     #else
1317     fmtopt[2] = 'd';
1318     #endif
1319 greg 2.1 if (sampcnt <= 0) sampcnt = 10000;
1320     }
1321     sprintf(sampcntbuf, "%d", sampcnt);
1322     CHECKARGC(3); /* add our format & sample count */
1323     rcarg[nrcargs++] = fmtopt;
1324     rcarg[nrcargs++] = "-c";
1325     rcarg[nrcargs++] = sampcntbuf;
1326     /* add receiver arguments to rcontrib */
1327     if (load_scene(argv[a], add_recv_object) < 0)
1328     return(1);
1329     finish_receiver();
1330     if (sendfn == NULL) { /* pass-through mode? */
1331     CHECKARGC(1); /* add octree */
1332     rcarg[nrcargs++] = oconv_command(argc-a, argv+a);
1333     rcarg[nrcargs] = NULL;
1334     return(my_exec(rcarg)); /* rcontrib does everything */
1335     }
1336     clear_params(&curparams, 0); /* else load sender surface & params */
1337 greg 2.24 curmod[0] = '\0';
1338 greg 2.1 if (load_scene(sendfn, add_send_object) < 0)
1339     return(1);
1340     if ((nsbins = prepare_sampler()) <= 0)
1341     return(1);
1342     CHECKARGC(3); /* add row count and octree */
1343     rcarg[nrcargs++] = "-y";
1344     sprintf(nsbinbuf, "%d", nsbins);
1345     rcarg[nrcargs++] = nsbinbuf;
1346     rcarg[nrcargs++] = oconv_command(argc-a, argv+a);
1347     rcarg[nrcargs] = NULL;
1348     /* open pipe to rcontrib process */
1349     if ((rcfp = popen_arglist(rcarg, "w")) == NULL)
1350     return(1);
1351     SET_FILE_BINARY(rcfp);
1352     #ifdef getc_unlocked
1353     flockfile(rcfp);
1354     #endif
1355 greg 2.22 if (verbose > 0) {
1356 greg 2.1 fprintf(stderr, "%s: sampling %d directions", progname, nsbins);
1357     if (curparams.nsurfs > 1)
1358 greg 2.7 fprintf(stderr, " (%d elements)\n", curparams.nsurfs);
1359 greg 2.1 else
1360     fputc('\n', stderr);
1361     }
1362     for (i = 0; i < nsbins; i++) /* send rcontrib ray samples */
1363     if (!(*curparams.sample_basis)(&curparams, i, rcfp))
1364     return(1);
1365 greg 2.23 return(pclose(rcfp) < 0); /* all finished! */
1366 greg 2.1 userr:
1367     if (a < argc-2)
1368     fprintf(stderr, "%s: unsupported option '%s'", progname, argv[a]);
1369 greg 2.13 fprintf(stderr, "Usage: %s [-v][rcontrib options] sender.rad receiver.rad [-i system.oct] [system.rad ..]\n",
1370 greg 2.1 progname);
1371     return(1);
1372     }