ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/util/rfluxmtx.c
Revision: 2.12
Committed: Mon Jul 28 16:30:18 2014 UTC (9 years, 9 months ago) by greg
Content type: text/plain
Branch: MAIN
CVS Tags: rad4R2
Changes since 2.11: +2 -1 lines
Log Message:
Added inclusion of rtprocess.h

File Contents

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