ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/util/rfluxmtx.c
Revision: 2.34
Committed: Tue Feb 2 18:02:32 2016 UTC (8 years, 2 months ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 2.33: +2 -58 lines
Log Message:
Moved declaration of popen to paths.h and put convert_command() into module

File Contents

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