ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/util/rfluxmtx.c
Revision: 2.50
Committed: Tue Dec 10 19:18:43 2019 UTC (4 years, 4 months ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 2.49: +2 -2 lines
Log Message: 
Fixed glitch in previous change

File Contents

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