ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/util/rfluxmtx.c
Revision: 2.19
Committed: Thu Dec 4 05:26:28 2014 UTC (9 years, 4 months ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 2.18: +3 -12 lines
Log Message:
Improved behavior of anisotropic reflections

File Contents

# Content
1 #ifndef lint
2 static const char RCSid[] = "$Id: rfluxmtx.c,v 2.18 2014/08/26 19:10:10 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 binv = "0"; /* uniform sampling -- one bin */
421 uniform = 1;
422 } else if (tolower(curparams.hemis[0]) == 's' &&
423 tolower(curparams.hemis[1]) == 'c') {
424 /* assign parameters */
425 if (curparams.hsiz <= 1) {
426 fputs(progname, stderr);
427 fputs(": missing size for Shirley-Chiu sampling!\n", stderr);
428 exit(1);
429 }
430 calfn = shirchiufn; shirchiufn = NULL;
431 sprintf(sbuf, "SCdim=%d,rNx=%g,rNy=%g,rNz=%g,Ux=%g,Uy=%g,Uz=%g",
432 curparams.hsiz,
433 curparams.nrm[0], curparams.nrm[1], curparams.nrm[2],
434 curparams.vup[0], curparams.vup[1], curparams.vup[2]);
435 params = savqstr(sbuf);
436 binv = "scbin";
437 nbins = "SCdim*SCdim";
438 } else if ((tolower(curparams.hemis[0]) == 'r') |
439 (tolower(curparams.hemis[0]) == 't')) {
440 calfn = reinhfn; reinhfn = NULL;
441 sprintf(sbuf, "MF=%d,rNx=%g,rNy=%g,rNz=%g,Ux=%g,Uy=%g,Uz=%g",
442 curparams.hsiz,
443 curparams.nrm[0], curparams.nrm[1], curparams.nrm[2],
444 curparams.vup[0], curparams.vup[1], curparams.vup[2]);
445 params = savqstr(sbuf);
446 binv = "rbin";
447 nbins = "Nrbins";
448 } else if (tolower(curparams.hemis[0]) == 'k' &&
449 !curparams.hemis[1] |
450 (tolower(curparams.hemis[1]) == 'f') |
451 (curparams.hemis[1] == '1')) {
452 calfn = kfullfn; kfullfn = NULL;
453 binf = "kbin";
454 nbins = "Nkbins";
455 } else if (tolower(curparams.hemis[0]) == 'k' &&
456 (tolower(curparams.hemis[1]) == 'h') |
457 (curparams.hemis[1] == '2')) {
458 calfn = khalffn; khalffn = NULL;
459 binf = "khbin";
460 nbins = "Nkhbins";
461 } else if (tolower(curparams.hemis[0]) == 'k' &&
462 (tolower(curparams.hemis[1]) == 'q') |
463 (curparams.hemis[1] == '4')) {
464 calfn = kquarterfn; kquarterfn = NULL;
465 binf = "kqbin";
466 nbins = "Nkqbins";
467 } else {
468 fprintf(stderr, "%s: unrecognized hemisphere sampling: h=%s\n",
469 progname, curparams.hemis);
470 exit(1);
471 }
472 if (!uniform & (curparams.slist->styp == ST_SOURCE)) {
473 SURF *sp;
474 for (sp = curparams.slist; sp != NULL; sp = sp->next)
475 if (fabs(sp->area - PI) > 1e-3) {
476 fprintf(stderr, "%s: source '%s' must be 180-degrees\n",
477 progname, sp->sname);
478 exit(1);
479 }
480 }
481 if (calfn != NULL) { /* add cal file if needed */
482 CHECKARGC(2);
483 rcarg[nrcargs++] = "-f";
484 rcarg[nrcargs++] = calfn;
485 }
486 if (params != NULL) { /* parameters _after_ cal file */
487 CHECKARGC(2);
488 rcarg[nrcargs++] = "-p";
489 rcarg[nrcargs++] = params;
490 }
491 if (nbins != NULL) { /* add #bins if set */
492 CHECKARGC(2);
493 rcarg[nrcargs++] = "-bn";
494 rcarg[nrcargs++] = nbins;
495 }
496 if (binv != NULL) {
497 CHECKARGC(2); /* assign bin variable */
498 rcarg[nrcargs++] = "-b";
499 rcarg[nrcargs++] = binv;
500 } else if (binf != NULL) {
501 CHECKARGC(2); /* assign bin function */
502 rcarg[nrcargs++] = "-b";
503 sprintf(sbuf, "%s(%g,%g,%g,%g,%g,%g)", binf,
504 curparams.nrm[0], curparams.nrm[1], curparams.nrm[2],
505 curparams.vup[0], curparams.vup[1], curparams.vup[2]);
506 rcarg[nrcargs++] = savqstr(sbuf);
507 }
508 CHECKARGC(2); /* modifier argument goes last */
509 rcarg[nrcargs++] = "-m";
510 rcarg[nrcargs++] = savqstr(curmod);
511 }
512
513 /* Make randomly oriented tangent plane axes for given normal direction */
514 static void
515 make_axes(FVECT uva[2], const FVECT nrm)
516 {
517 int i;
518
519 if (!getperpendicular(uva[0], nrm)) {
520 fputs(progname, stderr);
521 fputs(": bad surface normal in make_axes!\n", stderr);
522 exit(1);
523 }
524 fcross(uva[1], nrm, uva[0]);
525 }
526
527 /* Illegal sender surfaces end up here */
528 static int
529 ssamp_bad(FVECT orig, SURF *sp, double x)
530 {
531 fprintf(stderr, "%s: illegal sender surface '%s'\n",
532 progname, sp->sname);
533 return(0);
534 }
535
536 /* Generate origin on ring surface from uniform random variable */
537 static int
538 ssamp_ring(FVECT orig, SURF *sp, double x)
539 {
540 FVECT *uva = (FVECT *)sp->priv;
541 double samp2[2];
542 double uv[2];
543 int i;
544
545 if (uva == NULL) { /* need tangent axes */
546 uva = (FVECT *)malloc(sizeof(FVECT)*2);
547 if (uva == NULL) {
548 fputs(progname, stderr);
549 fputs(": out of memory in ssamp_ring!\n", stderr);
550 return(0);
551 }
552 make_axes(uva, sp->snrm);
553 sp->priv = (void *)uva;
554 }
555 SDmultiSamp(samp2, 2, x);
556 samp2[0] = sqrt(samp2[0]*sp->area*(1./PI) + sp->farg[6]*sp->farg[6]);
557 samp2[1] *= 2.*PI;
558 uv[0] = samp2[0]*tcos(samp2[1]);
559 uv[1] = samp2[0]*tsin(samp2[1]);
560 for (i = 3; i--; )
561 orig[i] = sp->farg[i] + uv[0]*uva[0][i] + uv[1]*uva[1][i];
562 return(1);
563 }
564
565 /* Add triangle to polygon's list (call-back function) */
566 static int
567 add_triangle(const Vert2_list *tp, int a, int b, int c)
568 {
569 POLYTRIS *ptp = (POLYTRIS *)tp->p;
570 struct ptri *trip = ptp->tri + ptp->ntris++;
571
572 trip->vndx[0] = a;
573 trip->vndx[1] = b;
574 trip->vndx[2] = c;
575 return(1);
576 }
577
578 /* Generate origin on polygon surface from uniform random variable */
579 static int
580 ssamp_poly(FVECT orig, SURF *sp, double x)
581 {
582 POLYTRIS *ptp = (POLYTRIS *)sp->priv;
583 double samp2[2];
584 double *v0, *v1, *v2;
585 int i;
586
587 if (ptp == NULL) { /* need to triangulate */
588 ptp = (POLYTRIS *)malloc(sizeof(POLYTRIS) +
589 sizeof(struct ptri)*(sp->nfargs/3 - 3));
590 if (ptp == NULL)
591 goto memerr;
592 if (sp->nfargs == 3) { /* simple case */
593 ptp->ntris = 1;
594 ptp->tri[0].vndx[0] = 0;
595 ptp->tri[0].vndx[1] = 1;
596 ptp->tri[0].vndx[2] = 2;
597 ptp->tri[0].afrac = 1;
598 } else {
599 Vert2_list *v2l = polyAlloc(sp->nfargs/3);
600 if (v2l == NULL)
601 goto memerr;
602 make_axes(ptp->uva, sp->snrm);
603 for (i = v2l->nv; i--; ) {
604 v2l->v[i].mX = DOT(sp->farg+3*i, ptp->uva[0]);
605 v2l->v[i].mY = DOT(sp->farg+3*i, ptp->uva[1]);
606 }
607 ptp->ntris = 0;
608 v2l->p = (void *)ptp;
609 if (!polyTriangulate(v2l, add_triangle)) {
610 fprintf(stderr,
611 "%s: cannot triangulate polygon '%s'\n",
612 progname, sp->sname);
613 return(0);
614 }
615 for (i = ptp->ntris; i--; ) {
616 int a = ptp->tri[i].vndx[0];
617 int b = ptp->tri[i].vndx[1];
618 int c = ptp->tri[i].vndx[2];
619 ptp->tri[i].afrac =
620 (v2l->v[a].mX*v2l->v[b].mY -
621 v2l->v[b].mX*v2l->v[a].mY +
622 v2l->v[b].mX*v2l->v[c].mY -
623 v2l->v[c].mX*v2l->v[b].mY +
624 v2l->v[c].mX*v2l->v[a].mY -
625 v2l->v[a].mX*v2l->v[c].mY) /
626 (2.*sp->area);
627 }
628 polyFree(v2l);
629 }
630 sp->priv = (void *)ptp;
631 }
632 /* pick triangle by partial area */
633 for (i = 0; i < ptp->ntris && x > ptp->tri[i].afrac; i++)
634 x -= ptp->tri[i].afrac;
635 SDmultiSamp(samp2, 2, x/ptp->tri[i].afrac);
636 samp2[0] *= samp2[1] = sqrt(samp2[1]);
637 samp2[1] = 1. - samp2[1];
638 v0 = sp->farg + 3*ptp->tri[i].vndx[0];
639 v1 = sp->farg + 3*ptp->tri[i].vndx[1];
640 v2 = sp->farg + 3*ptp->tri[i].vndx[2];
641 for (i = 3; i--; )
642 orig[i] = v0[i] + samp2[0]*(v1[i] - v0[i])
643 + samp2[1]*(v2[i] - v0[i]) ;
644 return(1);
645 memerr:
646 fputs(progname, stderr);
647 fputs(": out of memory in ssamp_poly!\n", stderr);
648 return(0);
649 }
650
651 /* Compute sample origin based on projected areas of sender subsurfaces */
652 static int
653 sample_origin(PARAMS *p, FVECT orig, const FVECT rdir, double x)
654 {
655 static double *projsa;
656 static int nall;
657 double tarea = 0;
658 int i;
659 SURF *sp;
660 /* special case for lone surface */
661 if (p->nsurfs == 1) {
662 sp = p->slist;
663 if (DOT(sp->snrm, rdir) >= -FTINY) {
664 fprintf(stderr,
665 "%s: internal - sample behind sender '%s'\n",
666 progname, sp->sname);
667 return(0);
668 }
669 return((*orig_in_surf[sp->styp])(orig, sp, x));
670 }
671 if (p->nsurfs > nall) { /* (re)allocate surface area cache */
672 if (projsa) free(projsa);
673 projsa = (double *)malloc(sizeof(double)*p->nsurfs);
674 if (!projsa) return(0);
675 nall = p->nsurfs;
676 }
677 /* compute projected areas */
678 for (i = 0, sp = p->slist; sp != NULL; i++, sp = sp->next) {
679 projsa[i] = -DOT(sp->snrm, rdir) * sp->area;
680 tarea += projsa[i] *= (double)(projsa[i] > FTINY);
681 }
682 if (tarea <= FTINY) { /* wrong side of sender? */
683 fputs(progname, stderr);
684 fputs(": internal - sample behind all sender elements!\n",
685 stderr);
686 return(0);
687 }
688 tarea *= x; /* get surface from list */
689 for (i = 0, sp = p->slist; tarea > projsa[i]; sp = sp->next)
690 tarea -= projsa[i++];
691 return((*orig_in_surf[sp->styp])(orig, sp, tarea/projsa[i]));
692 }
693
694 /* Uniform sample generator */
695 static int
696 sample_uniform(PARAMS *p, int b, FILE *fp)
697 {
698 int n = sampcnt;
699 double samp3[3];
700 double duvw[3];
701 FVECT orig_dir[2];
702 int i;
703
704 if (fp == NULL) /* just requesting number of bins? */
705 return(1);
706
707 while (n--) { /* stratified hemisphere sampling */
708 SDmultiSamp(samp3, 3, (n+frandom())/sampcnt);
709 SDsquare2disk(duvw, samp3[1], samp3[2]);
710 duvw[2] = -sqrt(1. - duvw[0]*duvw[0] - duvw[1]*duvw[1]);
711 for (i = 3; i--; )
712 orig_dir[1][i] = duvw[0]*p->udir[i] +
713 duvw[1]*p->vdir[i] +
714 duvw[2]*p->nrm[i] ;
715 if (!sample_origin(p, orig_dir[0], orig_dir[1], samp3[0]))
716 return(0);
717 if (fwrite(orig_dir, sizeof(FVECT), 2, fp) != 2)
718 return(0);
719 }
720 return(1);
721 }
722
723 /* Shirly-Chiu sample generator */
724 static int
725 sample_shirchiu(PARAMS *p, int b, FILE *fp)
726 {
727 int n = sampcnt;
728 double samp3[3];
729 double duvw[3];
730 FVECT orig_dir[2];
731 int i;
732
733 if (fp == NULL) /* just requesting number of bins? */
734 return(p->hsiz*p->hsiz);
735
736 while (n--) { /* stratified sampling */
737 SDmultiSamp(samp3, 3, (n+frandom())/sampcnt);
738 SDsquare2disk(duvw, (b/p->hsiz + samp3[1])/curparams.hsiz,
739 (b%p->hsiz + samp3[2])/curparams.hsiz);
740 duvw[2] = sqrt(1. - duvw[0]*duvw[0] - duvw[1]*duvw[1]);
741 for (i = 3; i--; )
742 orig_dir[1][i] = -duvw[0]*p->udir[i] -
743 duvw[1]*p->vdir[i] -
744 duvw[2]*p->nrm[i] ;
745 if (!sample_origin(p, orig_dir[0], orig_dir[1], samp3[0]))
746 return(0);
747 if (fwrite(orig_dir, sizeof(FVECT), 2, fp) != 2)
748 return(0);
749 }
750 return(1);
751 }
752
753 /* Reinhart/Tregenza sample generator */
754 static int
755 sample_reinhart(PARAMS *p, int b, FILE *fp)
756 {
757 #define T_NALT 7
758 static const int tnaz[T_NALT] = {30, 30, 24, 24, 18, 12, 6};
759 const int RowMax = T_NALT*p->hsiz + 1;
760 const double RAH = (.5*PI)/(RowMax-.5);
761 #define rnaz(r) (r >= RowMax-1 ? 1 : p->hsiz*tnaz[r/p->hsiz])
762 int n = sampcnt;
763 int row, col;
764 double samp3[3];
765 double alt, azi;
766 double duvw[3];
767 FVECT orig_dir[2];
768 int i;
769
770 if (fp == NULL) { /* just requesting number of bins? */
771 n = 0;
772 for (row = RowMax; row--; ) n += rnaz(row);
773 return(n);
774 }
775 row = 0; /* identify row & column */
776 col = b;
777 while (col >= rnaz(row)) {
778 col -= rnaz(row);
779 ++row;
780 }
781 while (n--) { /* stratified sampling */
782 SDmultiSamp(samp3, 3, (n+frandom())/sampcnt);
783 alt = (row+samp3[1])*RAH;
784 azi = (2.*PI)*(col+samp3[2]-.5)/rnaz(row);
785 duvw[2] = cos(alt); /* measured from horizon */
786 duvw[0] = tsin(azi)*duvw[2];
787 duvw[1] = tcos(azi)*duvw[2];
788 duvw[2] = sqrt(1. - duvw[2]*duvw[2]);
789 for (i = 3; i--; )
790 orig_dir[1][i] = -duvw[0]*p->udir[i] -
791 duvw[1]*p->vdir[i] -
792 duvw[2]*p->nrm[i] ;
793 if (!sample_origin(p, orig_dir[0], orig_dir[1], samp3[0]))
794 return(0);
795 if (fwrite(orig_dir, sizeof(FVECT), 2, fp) != 2)
796 return(0);
797 }
798 return(1);
799 #undef rnaz
800 #undef T_NALT
801 }
802
803 /* Klems sample generator */
804 static int
805 sample_klems(PARAMS *p, int b, FILE *fp)
806 {
807 static const char bname[4][20] = {
808 "LBNL/Klems Full",
809 "LBNL/Klems Half",
810 "INTERNAL ERROR",
811 "LBNL/Klems Quarter"
812 };
813 static ANGLE_BASIS *kbasis[4];
814 const int bi = p->hemis[1] - '1';
815 int n = sampcnt;
816 double samp2[2];
817 double duvw[3];
818 FVECT orig_dir[2];
819 int i;
820
821 if (!kbasis[bi]) { /* need to get basis, first */
822 for (i = 4; i--; )
823 if (!strcasecmp(abase_list[i].name, bname[bi])) {
824 kbasis[bi] = &abase_list[i];
825 break;
826 }
827 if (i < 0) {
828 fprintf(stderr, "%s: unknown hemisphere basis '%s'\n",
829 progname, bname[bi]);
830 return(0);
831 }
832 }
833 if (fp == NULL) /* just requesting number of bins? */
834 return(kbasis[bi]->nangles);
835
836 while (n--) { /* stratified sampling */
837 SDmultiSamp(samp2, 2, (n+frandom())/sampcnt);
838 if (!bi_getvec(duvw, b+samp2[1], kbasis[bi]))
839 return(0);
840 for (i = 3; i--; )
841 orig_dir[1][i] = duvw[0]*p->udir[i] +
842 duvw[1]*p->vdir[i] +
843 duvw[2]*p->nrm[i] ;
844 if (!sample_origin(p, orig_dir[0], orig_dir[1], samp2[0]))
845 return(0);
846 if (fwrite(orig_dir, sizeof(FVECT), 2, fp) != 2)
847 return(0);
848 }
849 return(1);
850 }
851
852 /* Prepare hemisphere basis sampler that will send rays to rcontrib */
853 static int
854 prepare_sampler(void)
855 {
856 if (curparams.slist == NULL) { /* missing sample surface! */
857 fputs(progname, stderr);
858 fputs(": no sender surface!\n", stderr);
859 return(-1);
860 }
861 /* misplaced output file spec. */
862 if ((curparams.outfn != NULL) & (verbose >= 0))
863 fprintf(stderr, "%s: warning - ignoring output file in sender ('%s')\n",
864 progname, curparams.outfn);
865 /* check/set basis hemisphere */
866 if (!curparams.hemis[0]) {
867 fputs(progname, stderr);
868 fputs(": missing sender sampling type!\n", stderr);
869 return(-1);
870 }
871 if (normalize(curparams.nrm) == 0) {
872 fputs(progname, stderr);
873 fputs(": undefined normal for sender sampling\n", stderr);
874 return(-1);
875 }
876 if (normalize(curparams.vup) == 0) {
877 if (fabs(curparams.nrm[2]) < .7)
878 curparams.vup[2] = 1;
879 else
880 curparams.vup[1] = 1;
881 }
882 VCROSS(curparams.udir, curparams.vup, curparams.nrm);
883 if (normalize(curparams.udir) == 0) {
884 fputs(progname, stderr);
885 fputs(": up vector coincides with sender normal\n", stderr);
886 return(-1);
887 }
888 VCROSS(curparams.vdir, curparams.nrm, curparams.udir);
889 if (tolower(curparams.hemis[0]) == 'u' | curparams.hemis[0] == '1')
890 curparams.sample_basis = sample_uniform;
891 else if (tolower(curparams.hemis[0]) == 's' &&
892 tolower(curparams.hemis[1]) == 'c')
893 curparams.sample_basis = sample_shirchiu;
894 else if ((tolower(curparams.hemis[0]) == 'r') |
895 (tolower(curparams.hemis[0]) == 't'))
896 curparams.sample_basis = sample_reinhart;
897 else if (tolower(curparams.hemis[0]) == 'k') {
898 switch (curparams.hemis[1]) {
899 case '1':
900 case '2':
901 case '4':
902 break;
903 case 'f':
904 case 'F':
905 case '\0':
906 curparams.hemis[1] = '1';
907 break;
908 case 'h':
909 case 'H':
910 curparams.hemis[1] = '2';
911 break;
912 case 'q':
913 case 'Q':
914 curparams.hemis[1] = '4';
915 break;
916 default:
917 goto unrecognized;
918 }
919 curparams.hemis[2] = '\0';
920 curparams.sample_basis = sample_klems;
921 } else
922 goto unrecognized;
923 /* return number of bins */
924 return((*curparams.sample_basis)(&curparams,0,NULL));
925 unrecognized:
926 fprintf(stderr, "%s: unrecognized sender sampling: h=%s\n",
927 progname, curparams.hemis);
928 return(-1);
929 }
930
931 /* Compute normal and area for polygon */
932 static int
933 finish_polygon(SURF *p)
934 {
935 const int nv = p->nfargs / 3;
936 FVECT e1, e2, vc;
937 int i;
938
939 memset(p->snrm, 0, sizeof(FVECT));
940 VSUB(e1, p->farg+3, p->farg);
941 for (i = 2; i < nv; i++) {
942 VSUB(e2, p->farg+3*i, p->farg);
943 VCROSS(vc, e1, e2);
944 p->snrm[0] += vc[0];
945 p->snrm[1] += vc[1];
946 p->snrm[2] += vc[2];
947 VCOPY(e1, e2);
948 }
949 p->area = normalize(p->snrm)*0.5;
950 return(p->area > FTINY);
951 }
952
953 /* Add a surface to our current parameters */
954 static void
955 add_surface(int st, const char *oname, FILE *fp)
956 {
957 SURF *snew;
958 int n;
959 /* get floating-point arguments */
960 if (!fscanf(fp, "%d", &n)) return;
961 while (n-- > 0) fscanf(fp, "%*s");
962 if (!fscanf(fp, "%d", &n)) return;
963 while (n-- > 0) fscanf(fp, "%*d");
964 if (!fscanf(fp, "%d", &n) || n <= 0) return;
965 snew = (SURF *)malloc(sizeof(SURF) + sizeof(double)*(n-1));
966 if (snew == NULL) {
967 fputs(progname, stderr);
968 fputs(": out of memory!\n", stderr);
969 exit(1);
970 }
971 strncpy(snew->sname, oname, sizeof(snew->sname)-1);
972 snew->sname[sizeof(snew->sname)-1] = '\0';
973 snew->styp = st;
974 snew->priv = NULL;
975 snew->nfargs = n;
976 for (n = 0; n < snew->nfargs; n++)
977 if (fscanf(fp, "%lf", &snew->farg[n]) != 1) {
978 fprintf(stderr, "%s: error reading arguments for '%s'\n",
979 progname, oname);
980 exit(1);
981 }
982 switch (st) {
983 case ST_RING:
984 if (snew->nfargs != 8)
985 goto badcount;
986 VCOPY(snew->snrm, snew->farg+3);
987 if (normalize(snew->snrm) == 0)
988 goto badnorm;
989 if (snew->farg[7] < snew->farg[6]) {
990 double t = snew->farg[7];
991 snew->farg[7] = snew->farg[6];
992 snew->farg[6] = t;
993 }
994 snew->area = PI*(snew->farg[7]*snew->farg[7] -
995 snew->farg[6]*snew->farg[6]);
996 break;
997 case ST_POLY:
998 if (snew->nfargs < 9 || snew->nfargs % 3)
999 goto badcount;
1000 finish_polygon(snew);
1001 break;
1002 case ST_SOURCE:
1003 if (snew->nfargs != 4)
1004 goto badcount;
1005 for (n = 3; n--; ) /* need to reverse "normal" */
1006 snew->snrm[n] = -snew->farg[n];
1007 if (normalize(snew->snrm) == 0)
1008 goto badnorm;
1009 snew->area = sin((PI/180./2.)*snew->farg[3]);
1010 snew->area *= PI*snew->area;
1011 break;
1012 }
1013 if ((snew->area <= FTINY) & (verbose >= 0)) {
1014 fprintf(stderr, "%s: warning - zero area for surface '%s'\n",
1015 progname, oname);
1016 free(snew);
1017 return;
1018 }
1019 VSUM(curparams.nrm, curparams.nrm, snew->snrm, snew->area);
1020 snew->next = curparams.slist;
1021 curparams.slist = snew;
1022 curparams.nsurfs++;
1023 return;
1024 badcount:
1025 fprintf(stderr, "%s: bad argument count for surface element '%s'\n",
1026 progname, oname);
1027 exit(1);
1028 badnorm:
1029 fprintf(stderr, "%s: bad orientation for surface element '%s'\n",
1030 progname, oname);
1031 exit(1);
1032 }
1033
1034 /* Parse a receiver object (look for modifiers to add to rcontrib command) */
1035 static int
1036 add_recv_object(FILE *fp)
1037 {
1038 int st;
1039 char thismod[128], otype[32], oname[128];
1040 int n;
1041
1042 if (fscanf(fp, "%s %s %s", thismod, otype, oname) != 3)
1043 return(0); /* must have hit EOF! */
1044 if (!strcmp(otype, "alias")) {
1045 fscanf(fp, "%*s"); /* skip alias */
1046 return(0);
1047 }
1048 /* is it a new receiver? */
1049 if ((st = surf_type(otype)) != ST_NONE) {
1050 if (curparams.slist != NULL && (st == ST_SOURCE) ^
1051 (curparams.slist->styp == ST_SOURCE)) {
1052 fputs(progname, stderr);
1053 fputs(": cannot mix source/non-source receivers!\n", stderr);
1054 return(-1);
1055 }
1056 if (strcmp(thismod, curmod)) {
1057 if (curmod[0]) { /* output last receiver? */
1058 finish_receiver();
1059 clear_params(&curparams, 1);
1060 }
1061 parse_params(&curparams, newparams);
1062 newparams[0] = '\0';
1063 strcpy(curmod, thismod);
1064 }
1065 add_surface(st, oname, fp); /* read & store surface */
1066 return(1);
1067 }
1068 /* else skip arguments */
1069 if (!fscanf(fp, "%d", &n)) return(0);
1070 while (n-- > 0) fscanf(fp, "%*s");
1071 if (!fscanf(fp, "%d", &n)) return(0);
1072 while (n-- > 0) fscanf(fp, "%*d");
1073 if (!fscanf(fp, "%d", &n)) return(0);
1074 while (n-- > 0) fscanf(fp, "%*f");
1075 return(0);
1076 }
1077
1078 /* Parse a sender object */
1079 static int
1080 add_send_object(FILE *fp)
1081 {
1082 int st;
1083 char otype[32], oname[128];
1084 int n;
1085
1086 if (fscanf(fp, "%*s %s %s", otype, oname) != 2)
1087 return(0); /* must have hit EOF! */
1088 if (!strcmp(otype, "alias")) {
1089 fscanf(fp, "%*s"); /* skip alias */
1090 return(0);
1091 }
1092 /* is it a new surface? */
1093 if ((st = surf_type(otype)) != ST_NONE) {
1094 if (st == ST_SOURCE) {
1095 fputs(progname, stderr);
1096 fputs(": cannot use source as a sender!\n", stderr);
1097 return(-1);
1098 }
1099 parse_params(&curparams, newparams);
1100 newparams[0] = '\0';
1101 add_surface(st, oname, fp); /* read & store surface */
1102 return(0);
1103 }
1104 /* else skip arguments */
1105 if (!fscanf(fp, "%d", &n)) return(0);
1106 while (n-- > 0) fscanf(fp, "%*s");
1107 if (!fscanf(fp, "%d", &n)) return(0);
1108 while (n-- > 0) fscanf(fp, "%*d");
1109 if (!fscanf(fp, "%d", &n)) return(0);
1110 while (n-- > 0) fscanf(fp, "%*f");
1111 return(0);
1112 }
1113
1114 /* Load a Radiance scene using the given callback function for objects */
1115 static int
1116 load_scene(const char *inspec, int (*ocb)(FILE *))
1117 {
1118 int rv = 0;
1119 char inpbuf[1024];
1120 FILE *fp;
1121 int c;
1122
1123 if (*inspec == '!')
1124 fp = popen(inspec+1, "r");
1125 else
1126 fp = fopen(inspec, "r");
1127 if (fp == NULL) {
1128 fprintf(stderr, "%s: cannot load '%s'\n", progname, inspec);
1129 return(-1);
1130 }
1131 while ((c = getc(fp)) != EOF) { /* load receiver data */
1132 if (isspace(c)) /* skip leading white space */
1133 continue;
1134 if (c == '!') { /* read from a new command */
1135 inpbuf[0] = c;
1136 if (fgetline(inpbuf+1, sizeof(inpbuf)-1, fp) != NULL) {
1137 if ((c = load_scene(inpbuf, ocb)) < 0)
1138 return(c);
1139 rv += c;
1140 }
1141 continue;
1142 }
1143 if (c == '#') { /* parameters/comment */
1144 if ((c = getc(fp)) == EOF || ungetc(c,fp) == EOF)
1145 break;
1146 if (!isspace(c) && fscanf(fp, "%s", inpbuf) == 1 &&
1147 !strcmp(inpbuf, PARAMSTART)) {
1148 if (fgets(inpbuf, sizeof(inpbuf), fp) != NULL)
1149 strcat(newparams, inpbuf);
1150 continue;
1151 }
1152 while ((c = getc(fp)) != EOF && c != '\n')
1153 ; /* else skipping comment */
1154 continue;
1155 }
1156 ungetc(c, fp); /* else check object for receiver */
1157 c = (*ocb)(fp);
1158 if (c < 0)
1159 return(c);
1160 rv += c;
1161 }
1162 /* close our input stream */
1163 c = (*inspec == '!') ? pclose(fp) : fclose(fp);
1164 if (c != 0) {
1165 fprintf(stderr, "%s: error loading '%s'\n", progname, inspec);
1166 return(-1);
1167 }
1168 return(rv);
1169 }
1170
1171 /* Get command arguments and run program */
1172 int
1173 main(int argc, char *argv[])
1174 {
1175 char fmtopt[6] = "-faa"; /* default output is ASCII */
1176 char *xrs=NULL, *yrs=NULL, *ldopt=NULL;
1177 char *iropt = NULL;
1178 char *sendfn;
1179 char sampcntbuf[32], nsbinbuf[32];
1180 FILE *rcfp;
1181 int nsbins;
1182 int a, i;
1183 /* screen rcontrib options */
1184 progname = argv[0];
1185 for (a = 1; a < argc-2; a++) {
1186 int na;
1187 /* check for argument expansion */
1188 while ((na = expandarg(&argc, &argv, a)) > 0)
1189 ;
1190 if (na < 0) {
1191 fprintf(stderr, "%s: cannot expand '%s'\n",
1192 progname, argv[a]);
1193 return(1);
1194 }
1195 if (argv[a][0] != '-' || !argv[a][1])
1196 break;
1197 na = 1;
1198 switch (argv[a][1]) { /* !! Keep consistent !! */
1199 case 'v': /* verbose mode */
1200 verbose = 1;
1201 na = 0;
1202 continue;
1203 case 'f': /* special case for -fo, -ff, etc. */
1204 switch (argv[a][2]) {
1205 case '\0': /* cal file */
1206 goto userr;
1207 case 'o': /* force output */
1208 goto userr;
1209 case 'a': /* output format */
1210 case 'f':
1211 case 'd':
1212 case 'c':
1213 if (!(fmtopt[3] = argv[a][3]))
1214 fmtopt[3] = argv[a][2];
1215 fmtopt[2] = argv[a][2];
1216 na = 0;
1217 continue; /* will pass later */
1218 default:
1219 goto userr;
1220 }
1221 break;
1222 case 'x': /* x-resolution */
1223 xrs = argv[++a];
1224 na = 0;
1225 continue;
1226 case 'y': /* y-resolution */
1227 yrs = argv[++a];
1228 na = 0;
1229 continue;
1230 case 'c': /* number of samples */
1231 sampcnt = atoi(argv[++a]);
1232 if (sampcnt <= 0)
1233 goto userr;
1234 na = 0; /* we re-add this later */
1235 continue;
1236 case 'I': /* only for pass-through mode */
1237 case 'i':
1238 iropt = argv[a];
1239 na = 0;
1240 continue;
1241 case 'w': /* options without arguments */
1242 if (argv[a][2] != '+') verbose = -1;
1243 case 'V':
1244 case 'u':
1245 case 'h':
1246 case 'r':
1247 break;
1248 case 'n': /* options with 1 argument */
1249 case 's':
1250 case 'o':
1251 na = 2;
1252 break;
1253 case 'b': /* special case */
1254 if (argv[a][2] != 'v') goto userr;
1255 break;
1256 case 'l': /* special case */
1257 if (argv[a][2] == 'd') {
1258 ldopt = argv[a];
1259 na = 0;
1260 continue;
1261 }
1262 na = 2;
1263 break;
1264 case 'd': /* special case */
1265 if (argv[a][2] != 'v') na = 2;
1266 break;
1267 case 'a': /* special case */
1268 na = (argv[a][2] == 'v') ? 4 : 2;
1269 break;
1270 case 'm': /* special case */
1271 if (!argv[a][2]) goto userr;
1272 na = (argv[a][2] == 'e') | (argv[a][2] == 'a') ? 4 : 2;
1273 break;
1274 default: /* anything else is verbotten */
1275 goto userr;
1276 }
1277 if (na <= 0) continue;
1278 CHECKARGC(na); /* pass on option */
1279 rcarg[nrcargs++] = argv[a];
1280 while (--na) /* + arguments if any */
1281 rcarg[nrcargs++] = argv[++a];
1282 }
1283 if (a > argc-2)
1284 goto userr; /* check at end of options */
1285 sendfn = argv[a++]; /* assign sender & receiver inputs */
1286 if (sendfn[0] == '-') { /* user wants pass-through mode? */
1287 if (sendfn[1]) goto userr;
1288 sendfn = NULL;
1289 if (iropt) {
1290 CHECKARGC(1);
1291 rcarg[nrcargs++] = iropt;
1292 }
1293 if (xrs) {
1294 CHECKARGC(2);
1295 rcarg[nrcargs++] = "-x";
1296 rcarg[nrcargs++] = xrs;
1297 }
1298 if (yrs) {
1299 CHECKARGC(2);
1300 rcarg[nrcargs++] = "-y";
1301 rcarg[nrcargs++] = yrs;
1302 }
1303 if (ldopt) {
1304 CHECKARGC(1);
1305 rcarg[nrcargs++] = ldopt;
1306 }
1307 if (sampcnt <= 0) sampcnt = 1;
1308 } else { /* else in sampling mode */
1309 if (iropt) {
1310 fputs(progname, stderr);
1311 fputs(": -i, -I supported for pass-through only\n", stderr);
1312 return(1);
1313 }
1314 fmtopt[2] = (sizeof(RREAL)==sizeof(double)) ? 'd' : 'f';
1315 if (sampcnt <= 0) sampcnt = 10000;
1316 }
1317 sprintf(sampcntbuf, "%d", sampcnt);
1318 CHECKARGC(3); /* add our format & sample count */
1319 rcarg[nrcargs++] = fmtopt;
1320 rcarg[nrcargs++] = "-c";
1321 rcarg[nrcargs++] = sampcntbuf;
1322 /* add receiver arguments to rcontrib */
1323 if (load_scene(argv[a], add_recv_object) < 0)
1324 return(1);
1325 finish_receiver();
1326 if (sendfn == NULL) { /* pass-through mode? */
1327 CHECKARGC(1); /* add octree */
1328 rcarg[nrcargs++] = oconv_command(argc-a, argv+a);
1329 rcarg[nrcargs] = NULL;
1330 return(my_exec(rcarg)); /* rcontrib does everything */
1331 }
1332 clear_params(&curparams, 0); /* else load sender surface & params */
1333 if (load_scene(sendfn, add_send_object) < 0)
1334 return(1);
1335 if ((nsbins = prepare_sampler()) <= 0)
1336 return(1);
1337 CHECKARGC(3); /* add row count and octree */
1338 rcarg[nrcargs++] = "-y";
1339 sprintf(nsbinbuf, "%d", nsbins);
1340 rcarg[nrcargs++] = nsbinbuf;
1341 rcarg[nrcargs++] = oconv_command(argc-a, argv+a);
1342 rcarg[nrcargs] = NULL;
1343 /* open pipe to rcontrib process */
1344 if ((rcfp = popen_arglist(rcarg, "w")) == NULL)
1345 return(1);
1346 SET_FILE_BINARY(rcfp);
1347 #ifdef getc_unlocked
1348 flockfile(rcfp);
1349 #endif
1350 if (verbose) {
1351 fprintf(stderr, "%s: sampling %d directions", progname, nsbins);
1352 if (curparams.nsurfs > 1)
1353 fprintf(stderr, " (%d elements)\n", curparams.nsurfs);
1354 else
1355 fputc('\n', stderr);
1356 }
1357 for (i = 0; i < nsbins; i++) /* send rcontrib ray samples */
1358 if (!(*curparams.sample_basis)(&curparams, i, rcfp))
1359 return(1);
1360 return(pclose(rcfp) == 0); /* all finished! */
1361 userr:
1362 if (a < argc-2)
1363 fprintf(stderr, "%s: unsupported option '%s'", progname, argv[a]);
1364 fprintf(stderr, "Usage: %s [-v][rcontrib options] sender.rad receiver.rad [-i system.oct] [system.rad ..]\n",
1365 progname);
1366 return(1);
1367 }