ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/util/rfluxmtx.c
Revision: 2.47
Committed: Sun Mar 3 17:01:13 2019 UTC (5 years, 1 month ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 2.46: +2 -2 lines
Log Message:
Changed test to strict zero because it can sometimes be that!

File Contents

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