ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/util/rfluxmtx.c
Revision: 2.38
Committed: Thu Aug 18 00:52:48 2016 UTC (7 years, 7 months ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 2.37: +5 -5 lines
Log Message:
Switched over to more efficient fread/fwrite replacements getbinary/putbinary

File Contents

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