ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/util/rfluxmtx.c
Revision: 2.20
Committed: Mon Feb 2 16:04:11 2015 UTC (9 years, 1 month ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 2.19: +5 -2 lines
Log Message:
Fixed uniform receivers so backside does not get contributions

File Contents

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