ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/util/rcomb.c
Revision: 2.17
Committed: Thu May 23 17:13:52 2024 UTC (11 months, 1 week ago) by greg
Content type: text/plain
Branch: MAIN
Changes since 2.16: +8 -7 lines
Log Message:
fix(rcomb): Error in record size calc and cleanup of child processes on failure

File Contents

# User Rev Content
1 greg 2.1 #ifndef lint
2 greg 2.17 static const char RCSid[] = "$Id: rcomb.c,v 2.16 2024/05/23 15:48:44 greg Exp $";
3 greg 2.1 #endif
4     /*
5     * General component matrix combiner, operating on a row at a time.
6     */
7    
8     #include <math.h>
9     #include "platform.h"
10 greg 2.10 #include "rtprocess.h"
11 greg 2.1 #include "rtio.h"
12     #include "rmatrix.h"
13     #include "calcomp.h"
14    
15     #ifndef M_PI
16     #define M_PI 3.14159265358979323846
17     #endif
18    
19     #define MAXCOMP MAXCSAMP /* #components we support */
20    
21     /* Unary matrix operation(s) */
22     typedef struct {
23     double cmat[MAXCOMP*MAXCOMP]; /* component transformation */
24     double sca[MAXCOMP]; /* scalar coefficients */
25     const char *csym; /* symbolic coefficients */
26     short clen; /* number of coefficients */
27     short nsf; /* number of scalars */
28     } RUNARYOP;
29    
30     /* Input matrix */
31     typedef struct {
32     const char *inspec; /* input specification */
33     RUNARYOP preop; /* transform operation */
34     RMATRIX imx; /* input matrix header info */
35     RMATRIX *rmp; /* active single-row matrix */
36     FILE *infp; /* open input stream */
37     } ROPMAT;
38    
39     ROPMAT *mop = NULL; /* allocated input array */
40     int nall = 0; /* number allocated */
41     int nmats = 0; /* number of actual inputs */
42    
43     RMATRIX *mcat = NULL; /* final concatenation */
44     int mcat_last = 0; /* goes after trailing ops? */
45    
46 greg 2.2 int in_nrows; /* number of input rows (or 0) */
47     #define in_ncols (mop[0].rmp->ncols) /* number of input columns */
48 greg 2.1 #define in_ncomp (mop[0].rmp->ncomp) /* input #components */
49    
50     extern int nowarn; /* turn off warnings? */
51    
52     int cur_row; /* current input/output row */
53     int cur_col; /* current input/output column */
54     int cur_chan; /* if we're looping channels */
55    
56 greg 2.10 SUBPROC *cproc = NULL; /* child process array */
57     int nchildren = 0; /* # of child processes */
58     int inchild = -1; /* our child ID (-1: parent) */
59    
60 greg 2.1 static int checksymbolic(ROPMAT *rop);
61    
62     static int
63     split_input(ROPMAT *rop)
64     {
65     if (rop->rmp == &rop->imx && !(rop->rmp = rmx_copy(&rop->imx))) {
66     fputs("Out of memory in split_input()\n", stderr);
67     return(0);
68     }
69     rmx_reset(rop->rmp);
70     return(1);
71     }
72    
73     /* Check/set transform based on a reference input file */
74     static int
75     checkreffile(ROPMAT *rop)
76     {
77     static const char *curRF = NULL;
78     static RMATRIX refm;
79     const int nc = rop->imx.ncomp;
80     int i;
81    
82     if (!curRF || strcmp(rop->preop.csym, curRF)) {
83     FILE *fp = fopen(rop->preop.csym, "rb");
84     if (!rmx_load_header(&refm, fp)) {
85     fprintf(stderr, "%s: cannot read info header\n",
86     rop->preop.csym);
87     curRF = NULL;
88     if (fp) fclose(fp);
89     return(0);
90     }
91     fclose(fp);
92     curRF = rop->preop.csym;
93     }
94     if (refm.ncomp == 3) {
95     rop->preop.csym = (refm.dtype == DTxyze) ? "XYZ" : "RGB";
96     return(checksymbolic(rop));
97     }
98     if (refm.ncomp == 2) {
99     fprintf(stderr, "%s: cannot convert to 2 components\n",
100     curRF);
101     return(0);
102     }
103     if (refm.ncomp == 1) {
104     rop->preop.csym = "Y"; /* XXX big assumption */
105     return(checksymbolic(rop));
106     }
107     if (refm.ncomp == nc &&
108     !memcmp(refm.wlpart, rop->imx.wlpart, sizeof(refm.wlpart)))
109     return(1); /* nothing to do */
110    
111     if ((nc <= 3) | (nc > MAXCSAMP) | (refm.ncomp > MAXCSAMP)) {
112     fprintf(stderr, "%s: cannot resample from %d to %d components\n",
113     curRF, nc, refm.ncomp);
114     return(0);
115     }
116     if (!split_input(rop)) /* get our own struct */
117     return(0);
118     rop->preop.clen = refm.ncomp * nc; /* compute spec to ref */
119    
120     for (i = 0; i < nc; i++) {
121     SCOLOR scstim, scresp;
122     int j;
123     memset(scstim, 0, sizeof(COLORV)*nc);
124     scstim[i] = 1.f;
125     convertscolor(scresp, refm.ncomp, refm.wlpart[0], refm.wlpart[3],
126     scstim, nc, rop->imx.wlpart[0], rop->imx.wlpart[3]);
127     for (j = refm.ncomp; j-- > 0; )
128     rop->preop.cmat[j*nc + i] = scresp[j];
129     }
130     /* remember new spectral params */
131     memcpy(rop->rmp->wlpart, refm.wlpart, sizeof(rop->rmp->wlpart));
132     rop->rmp->ncomp = refm.ncomp;
133     return(1);
134     }
135    
136     /* Compute conversion row from spectrum to one channel of RGB */
137     static void
138     rgbrow(ROPMAT *rop, int r, int p)
139     {
140     const int nc = rop->imx.ncomp;
141     const float * wlp = rop->imx.wlpart;
142     int i;
143    
144     for (i = nc; i--; ) {
145     int nmEnd = wlp[0] + (wlp[3] - wlp[0])*i/nc;
146     int nmStart = wlp[0] + (wlp[3] - wlp[0])*(i+1)/nc;
147     COLOR crgb;
148     spec_rgb(crgb, nmStart, nmEnd);
149     rop->preop.cmat[r*nc+i] = crgb[p];
150     }
151     }
152    
153     /* Compute conversion row from spectrum to one channel of XYZ */
154     static void
155     xyzrow(ROPMAT *rop, int r, int p)
156     {
157     const int nc = rop->imx.ncomp;
158     const float * wlp = rop->imx.wlpart;
159     int i;
160    
161     for (i = nc; i--; ) {
162     int nmEnd = wlp[0] + (wlp[3] - wlp[0])*i/nc;
163     int nmStart = wlp[0] + (wlp[3] - wlp[0])*(i+1)/nc;
164     COLOR cxyz;
165     spec_cie(cxyz, nmStart, nmEnd);
166     rop->preop.cmat[r*nc+i] = cxyz[p];
167     }
168     }
169    
170     /* Use the spectral sensitivity function to compute matrix coefficients */
171     static void
172     sensrow(ROPMAT *rop, int r, double (*sf)(SCOLOR sc, int ncs, const float wlpt[4]))
173     {
174     const int nc = rop->imx.ncomp;
175     int i;
176    
177     for (i = nc; i--; ) {
178     SCOLOR sclr;
179     memset(sclr, 0, sizeof(COLORV)*nc);
180     sclr[i] = 1.f;
181     rop->preop.cmat[r*nc+i] = (*sf)(sclr, nc, rop->imx.wlpart);
182     }
183     }
184    
185     /* Check/set symbolic transform */
186     static int
187     checksymbolic(ROPMAT *rop)
188     {
189     const int nc = rop->imx.ncomp;
190     const int dt = rop->imx.dtype;
191 greg 2.3 double cf = 1;
192 greg 2.1 int i, j;
193     /* check suffix => reference file */
194     if (strchr(rop->preop.csym, '.') > rop->preop.csym)
195     return(checkreffile(rop));
196    
197     if (nc < 3) {
198     fprintf(stderr, "%s: -c '%s' requires at least 3 components\n",
199     rop->inspec, rop->preop.csym);
200     return(0);
201     }
202     rop->preop.clen = strlen(rop->preop.csym) * nc;
203     if (rop->preop.clen > MAXCOMP*MAXCOMP) {
204     fprintf(stderr, "%s: -c '%s' results in too many components\n",
205     rop->inspec, rop->preop.csym);
206     return(0);
207     }
208     for (j = 0; rop->preop.csym[j]; j++) {
209     int comp = 0;
210     switch (rop->preop.csym[j]) {
211     case 'B':
212 greg 2.3 case 'b':
213 greg 2.1 ++comp;
214     /* fall through */
215     case 'G':
216 greg 2.3 case 'g':
217 greg 2.1 ++comp;
218     /* fall through */
219     case 'R':
220 greg 2.3 case 'r':
221     if (rop->preop.csym[j] <= 'Z')
222     cf = 1./WHTEFFICACY;
223 greg 2.1 if (dt == DTxyze) {
224     for (i = 3; i--; )
225 greg 2.3 rop->preop.cmat[j*nc+i] = cf*xyz2rgbmat[comp][i];
226 greg 2.1 } else if (nc == 3)
227     rop->preop.cmat[j*nc+comp] = 1.;
228     else
229     rgbrow(rop, j, comp);
230     break;
231     case 'Z':
232 greg 2.3 case 'z':
233 greg 2.1 ++comp;
234     /* fall through */
235     case 'Y':
236 greg 2.3 case 'y':
237 greg 2.1 ++comp;
238     /* fall through */
239     case 'X':
240 greg 2.3 case 'x':
241     if ((rop->preop.csym[j] <= 'Z') & (dt != DTxyze))
242     cf = WHTEFFICACY;
243 greg 2.1 if (dt == DTxyze) {
244     rop->preop.cmat[j*nc+comp] = 1.;
245     } else if (nc == 3) {
246     for (i = 3; i--; )
247     rop->preop.cmat[j*nc+i] =
248     rgb2xyzmat[comp][i];
249     } else if (comp == CIEY)
250     sensrow(rop, j, scolor2photopic);
251     else
252     xyzrow(rop, j, comp);
253    
254 greg 2.3 for (i = nc*(cf != 1); i--; )
255     rop->preop.cmat[j*nc+i] *= cf;
256 greg 2.1 break;
257     case 'S': /* scotopic (il)luminance */
258 greg 2.3 cf = WHTSCOTOPIC;
259     /* fall through */
260     case 's':
261 greg 2.1 sensrow(rop, j, scolor2scotopic);
262 greg 2.3 for (i = nc*(cf != 1); i--; )
263     rop->preop.cmat[j*nc+i] *= cf;
264 greg 2.1 break;
265     case 'M': /* melanopic (il)luminance */
266 greg 2.3 cf = WHTMELANOPIC;
267     /* fall through */
268     case 'm':
269 greg 2.1 sensrow(rop, j, scolor2melanopic);
270 greg 2.3 for (i = nc*(cf != 1); i--; )
271     rop->preop.cmat[j*nc+i] *= cf;
272 greg 2.1 break;
273     case 'A': /* average component */
274 greg 2.3 case 'a':
275 greg 2.1 for (i = nc; i--; )
276     rop->preop.cmat[j*nc+i] = 1./(double)nc;
277     break;
278     default:
279     fprintf(stderr, "%s: -c '%c' unsupported\n",
280     rop->inspec, rop->preop.csym[j]);
281     return(0);
282     }
283     }
284     if (!split_input(rop)) /* get our own struct */
285     return(0);
286     memcpy(rop->rmp->wlpart, WLPART, sizeof(rop->rmp->wlpart));
287     rop->rmp->ncomp = rop->preop.clen / nc;
288     /* decide on output type */
289 greg 2.4 if (!strcasecmp(rop->preop.csym, "XYZ")) {
290 greg 2.1 if (dt <= DTspec)
291     rop->rmp->dtype = DTxyze;
292 greg 2.4 } else if (!strcasecmp(rop->preop.csym, "RGB")) {
293 greg 2.1 if (dt <= DTspec)
294     rop->rmp->dtype = DTrgbe;
295     } else if (rop->rmp->dtype == DTspec)
296     rop->rmp->dtype = DTfloat;
297     return(1);
298     }
299    
300     static int
301     get_component_xfm(ROPMAT *rop)
302     {
303     int i, j;
304    
305     if (rop->rmp != &rop->imx) { /* reset destination matrix */
306     rmx_free(rop->rmp);
307     rop->rmp = &rop->imx;
308     }
309     if (rop->preop.csym && /* symbolic transform? */
310     !checksymbolic(rop))
311     return(0);
312     /* undo exposure? */
313     if (fabs(1. - bright(rop->rmp->cexp)) > .025) {
314     if (rop->rmp->ncomp == 1)
315     rop->rmp->cexp[RED] = rop->rmp->cexp[GRN] =
316     rop->rmp->cexp[BLU] = bright(rop->rmp->cexp);
317     if (rop->preop.nsf <= 0) {
318     rop->preop.nsf = i = rop->rmp->ncomp;
319     while (i--)
320     rop->preop.sca[i] = 1.;
321     }
322     if (rop->preop.nsf == 1) {
323     if (rop->rmp->ncomp == 3) {
324     rop->preop.sca[2] = rop->preop.sca[1] =
325     rop->preop.sca[0];
326     rop->preop.nsf = 3;
327     } else
328     rop->preop.sca[0] /= bright(rop->rmp->cexp);
329     }
330     if (rop->preop.nsf == 3) {
331     opcolor(rop->preop.sca, /=, rop->rmp->cexp);
332     } else if (rop->preop.nsf > 3) { /* punt */
333     double mult = 1./bright(rop->rmp->cexp);
334     for (i = rop->preop.nsf; i--; )
335     rop->preop.sca[i] *= mult;
336     }
337     setcolor(rop->rmp->cexp, 1., 1., 1.);
338     }
339     if (rop->preop.clen > 0) { /* use component transform? */
340     if (rop->preop.clen % rop->imx.ncomp) {
341     fprintf(stderr, "%s: -c must have N x %d coefficients\n",
342     rop->inspec, rop->imx.ncomp);
343     return(0);
344     }
345     if (rop->preop.nsf > 0) { /* scale transform, instead */
346     if (rop->preop.nsf == 1) {
347     for (i = rop->preop.clen; i--; )
348     rop->preop.cmat[i] *= rop->preop.sca[0];
349     } else if (rop->preop.nsf*rop->imx.ncomp != rop->preop.clen) {
350     fprintf(stderr, "%s: -s must have one or %d factors\n",
351     rop->inspec,
352     rop->preop.clen/rop->imx.ncomp);
353     return(0);
354     } else {
355     for (i = rop->preop.nsf; i--; )
356     for (j = rop->imx.ncomp; j--; )
357     rop->preop.cmat[i*rop->imx.ncomp+j]
358     *= rop->preop.sca[i];
359     }
360     }
361     rop->preop.nsf = 0; /* now folded in */
362     if (!split_input(rop)) /* get our own struct */
363     return(0);
364     rop->rmp->ncomp = rop->preop.clen / rop->imx.ncomp;
365     if ((rop->rmp->ncomp > 3) & (rop->rmp->dtype <= DTspec)) {
366     rop->rmp->dtype = DTfloat; /* probably not actual spectrum */
367     memcpy(rop->rmp->wlpart, WLPART, sizeof(rop->rmp->wlpart));
368     }
369     } else if (rop->preop.nsf > 0) { /* else use scalar(s)? */
370     if (rop->preop.nsf == 1) {
371     for (i = rop->rmp->ncomp; --i; )
372     rop->preop.sca[i] = rop->preop.sca[0];
373     rop->preop.nsf = rop->rmp->ncomp;
374     } else if (rop->preop.nsf != rop->rmp->ncomp) {
375     fprintf(stderr, "%s: -s must have one or %d factors\n",
376     rop->inspec, rop->rmp->ncomp);
377     return(0);
378     }
379     }
380     return(1);
381     }
382    
383     static int
384     apply_op(RMATRIX *dst, const RMATRIX *src, const RUNARYOP *ro)
385     {
386     if (ro->clen > 0) {
387     RMATRIX *res = rmx_transform(src, dst->ncomp, ro->cmat);
388     if (!res) {
389     fputs("Error in call to rmx_transform()\n", stderr);
390     return(0);
391     }
392     if (!rmx_transfer_data(dst, res, 0))
393     return(0);
394     rmx_free(res);
395     } else if (dst != src)
396 greg 2.9 memcpy(dst->mtx, src->mtx, rmx_array_size(dst));
397 greg 2.1 if (ro->nsf == dst->ncomp)
398     rmx_scale(dst, ro->sca);
399     return(1);
400     }
401    
402     static int
403     open_input(ROPMAT *rop)
404     {
405     int outtype;
406    
407     if (!rop || !rop->inspec || !rop->inspec[0])
408     return(0);
409     if (rop->inspec == stdin_name)
410     rop->infp = stdin;
411     else if (rop->inspec[0] == '!')
412     rop->infp = popen(rop->inspec+1, "r");
413     else
414     rop->infp = fopen(rop->inspec, "rb");
415    
416     if (!rmx_load_header(&rop->imx, rop->infp)) {
417     fprintf(stderr, "Bad header from: %s\n", rop->inspec);
418     return(0);
419     }
420     return(get_component_xfm(rop));
421     }
422    
423     /* Return nominal wavelength associated with input component (return nm) */
424     static double
425     l_wavelength(char *nam)
426     {
427     double comp = argument(1);
428    
429     if ((comp < -.5) | (comp >= in_ncomp+.5)) {
430     errno = EDOM;
431     return(.0);
432     }
433     if (comp < .5) /* asking for #components? */
434     return(in_ncomp);
435    
436     if (in_ncomp == 3) { /* special case for RGB */
437     const int w0 = (int)(comp - .5);
438     return(mop[0].rmp->wlpart[w0] +
439     (comp-.5)*(mop[0].rmp->wlpart[w0+1] -
440     mop[0].rmp->wlpart[w0]));
441     }
442     return(mop[0].rmp->wlpart[0] + /* general case, even div. */
443     (comp-.5)/(double)in_ncomp *
444     (mop[0].rmp->wlpart[3] - mop[0].rmp->wlpart[0]));
445     }
446    
447     /* Return ith input with optional channel selector */
448     static double
449     l_chanin(char *nam)
450     {
451     double inp = argument(1);
452     int mi, chan;
453    
454     if ((mi = (int)(inp-.5)) < 0 || mi >= nmats) {
455     errno = EDOM;
456     return(.0);
457     }
458     if (inp < .5) /* asking for #inputs? */
459     return(nmats);
460    
461     if (nargum() >= 2) {
462     double cval = argument(2);
463     if (cval < .5 || (chan = (int)(cval-.5)) >= in_ncomp) {
464     errno = EDOM;
465     return(.0);
466     }
467     } else
468     chan = cur_chan;
469    
470     return(mop[mi].rmp->mtx[cur_col*in_ncomp + chan]);
471     }
472    
473     static int
474     initialize(RMATRIX *imp)
475     {
476     int i;
477     /* XXX struct is zeroed coming in */
478     setcolor(imp->cexp, 1.f, 1.f, 1.f);
479     for (i = 0; i < nmats; i++) { /* open each input */
480     int restype;
481     if (!open_input(&mop[i]))
482     return(0);
483     restype = mop[i].rmp->dtype;
484     if (!imp->dtype || (restype = rmx_newtype(restype, imp->dtype)) > 0)
485     imp->dtype = restype;
486     else
487     fprintf(stderr, "%s: warning - data type mismatch\n",
488     mop[i].inspec);
489     if (!i) {
490     imp->ncols = mop[0].rmp->ncols;
491     imp->ncomp = mop[0].rmp->ncomp;
492     memcpy(imp->wlpart, mop[0].rmp->wlpart, sizeof(imp->wlpart));
493 greg 2.2 } else if ((mop[i].rmp->ncols != imp->ncols) |
494     (mop[i].rmp->ncomp != imp->ncomp) |
495     ((in_nrows > 0) & (mop[i].rmp->nrows > 0) &
496     (mop[i].rmp->nrows != in_nrows))) {
497 greg 2.1 fprintf(stderr, "%s: mismatch in size or #components\n",
498     mop[i].inspec);
499     return(0);
500     } /* XXX should check wlpart? */
501 greg 2.2 if (in_nrows <= 0)
502     in_nrows = imp->nrows = mop[i].rmp->nrows;
503 greg 2.1 } /* set up .cal environment */
504     esupport |= E_VARIABLE|E_FUNCTION|E_RCONST;
505     esupport &= ~(E_OUTCHAN|E_INCHAN);
506     varset("PI", ':', M_PI);
507     varset("nfiles", ':', nmats);
508     varset("nrows", ':', in_nrows);
509     varset("ncols", ':', in_ncols);
510     varset("ncomp", ':', in_ncomp);
511     varset("R", ':', 1.);
512     varset("G", ':', 2.);
513     varset("B", ':', 3.);
514     funset("wl", 1, ':', l_wavelength);
515     funset("ci", 1, '=', l_chanin);
516     scompile("ri(i)=ci(i,R);gi(i)=ci(i,G);bi(i)=ci(i,B)", NULL, 0);
517     return(1);
518     }
519    
520     static void
521     output_headinfo(FILE *fp)
522     {
523     int i;
524    
525     for (i = 0; i < nmats; i++) {
526     const char *cp = mop[i].imx.info;
527     fputs(mop[i].inspec, fp);
528     fputs(":\n", fp);
529     if (!cp) continue;
530     while (*cp) {
531     if (*cp == '\n') {
532     cp++; /* avoid inadvertant terminus */
533     continue;
534     }
535     fputc('\t', fp); /* indent this input's info */
536     do
537     putc(*cp, fp);
538     while (*cp++ != '\n');
539     }
540     }
541     }
542    
543     static int
544 greg 2.10 spawned_children(int np)
545 greg 2.1 {
546 greg 2.17 size_t recsize = 0;
547 greg 2.10 int i, rv;
548    
549     #if defined(_WIN32) || defined(_WIN64)
550     if (np > 1) {
551     fputs("Warning: only one process under Windows\n", stderr);
552     np = 1;
553     } else
554     #endif
555     if ((in_nrows > 0) & (np > in_nrows))
556     np = in_nrows;
557     /* we'll be doing a row at a time */
558 greg 2.1 for (i = 0; i < nmats; i++) {
559 greg 2.10 mop[i].imx.nrows = 1;
560 greg 2.1 if (!rmx_prepare(&mop[i].imx))
561     goto memerror;
562 greg 2.17 recsize += rmx_array_size(&mop[i].imx);
563 greg 2.1 if (mop[i].rmp != &mop[i].imx) {
564     mop[i].rmp->nrows = 1;
565     if (!rmx_prepare(mop[i].rmp))
566     goto memerror;
567     }
568     }
569 greg 2.10 /* prep output row buffer */
570     if (mcat || mop[nmats].preop.clen > 0) {
571     if (!split_input(&mop[nmats])) /* need separate buffer */
572     return(0);
573     if (mop[nmats].preop.clen > 0)
574     mop[nmats].rmp->ncomp = mop[nmats].preop.clen /
575     mop[nmats].imx.ncomp;
576     mop[nmats].rmp->nrows = 1;
577     if (!mcat | !mcat_last && !rmx_prepare(mop[nmats].rmp))
578 greg 2.1 goto memerror;
579     }
580 greg 2.10 mop[nmats].imx.nrows = 1;
581     if (!rmx_prepare(&mop[nmats].imx))
582     goto memerror;
583     if (np <= 1) { /* single process return point */
584     #ifdef getc_unlocked
585     for (i = 0; i < nmats; i++)
586     flockfile(mop[i].infp);
587     flockfile(stdout);
588     #endif
589     return(0);
590     }
591     fflush(stdout); /* flush header & spawn children */
592     cproc = (SUBPROC *)malloc(sizeof(SUBPROC)*np);
593     if (!cproc)
594     goto memerror;
595     nchildren = np;
596     for (i = 0; i < np; i++) {
597     cproc[i].flags = PF_FILT_OUT;
598     cproc[i].w = dup(1);
599     cproc[i].r = 0;
600     cproc[i].pid = -1;
601     rv = open_process(&cproc[i], NULL);
602     if (rv <= 0) break;
603 greg 2.17 if (!i && 2*rv >= recsize) {
604 greg 2.16 fputs("Problem too small for multi-processing\n",
605     stderr);
606 greg 2.17 close_processes(cproc, 1);
607 greg 2.16 exit(1);
608     }
609 greg 2.10 }
610     if (rv < 0) {
611     perror("fork");
612 greg 2.17 close_processes(cproc, i);
613 greg 2.10 exit(1);
614     }
615 greg 2.17 if (rv > 0) /* parent return? */
616     return(1);
617 greg 2.11 inchild = i; /* our child index */
618     while (i-- > 0) /* don't share siblings' pipes */
619     close(cproc[i].w);
620 greg 2.10 fpurge(stdin); /* discard previous matrix input */
621     #ifdef getc_unlocked
622     flockfile(stdin);
623     #endif
624     for (i = 0; i < nmats; i++) {
625     if (mop[i].infp != stdin)
626 greg 2.11 fclose(mop[i].infp); /* ! pclose() */
627 greg 2.10 mop[i].infp = stdin;
628     mop[i].imx.dtype = DTdouble;
629     }
630     return(0); /* child return */
631     memerror:
632     fputs("Out of memory in spawned_children()\n", stderr);
633     exit(1);
634     }
635    
636     static int
637     parent_loop()
638     {
639     FILE **outfp = (FILE **)malloc(nchildren*sizeof(FILE *));
640     int i;
641    
642     if (!outfp) goto memerror;
643     for (i = 0; i < nchildren; i++) {
644     outfp[i] = fdopen(cproc[i].w, "w");
645     if (!outfp[i]) goto memerror;
646     #ifdef getc_unlocked
647     flockfile(outfp[i]);
648     #endif
649     }
650     #ifdef getc_unlocked
651     for (i = 0; i < nmats; i++)
652     flockfile(mop[i].infp);
653     #endif
654     for (cur_row = 0; (in_nrows <= 0) | (cur_row < in_nrows); cur_row++) {
655     FILE *ofp = outfp[cur_row % nchildren];
656 greg 2.11 for (i = 0; i < nmats; i++)
657 greg 2.10 if (!rmx_load_row(mop[i].imx.mtx, &mop[i].imx, mop[i].infp)) {
658     if (cur_row > in_nrows) /* unknown #input rows? */
659     break;
660     fprintf(stderr, "%s: read error at row %d\n",
661     mop[i].inspec, cur_row);
662     return(0);
663     }
664 greg 2.11 if (i < nmats)
665     break;
666     for (i = 0; i < nmats; i++)
667 greg 2.10 if (!rmx_write_data(mop[i].imx.mtx, mop[i].imx.ncomp,
668 greg 2.12 mop[i].imx.ncols, DTdouble, ofp))
669 greg 2.10 return(0);
670     if (fflush(ofp) == EOF)
671     return(0);
672     }
673 greg 2.16 for (i = 0; i < nchildren; i++) {
674     sleep(2); /* try to maintain order */
675 greg 2.10 fclose(outfp[i]);
676 greg 2.16 }
677 greg 2.10 free(outfp);
678     i = close_processes(cproc, nchildren);
679     free(cproc); cproc = NULL;
680     if (i < 0) {
681     fputs("Warning: missing child in parent_loop()\n", stderr);
682     return(1);
683     }
684     if (i > 0) {
685     fprintf(stderr, "Child exited with status %d\n", i);
686     return(0);
687     }
688     return(1);
689     memerror:
690     fputs("Out of memory in parent_loop()\n", stderr);
691     exit(1);
692     }
693    
694     static int
695     combine_input()
696     {
697     const int row0 = (inchild >= 0)*inchild;
698     const int rstep = nchildren + !nchildren;
699     ROPMAT *res = &mop[nmats];
700     int set_r, set_c;
701     RMATRIX *tmp = NULL;
702     int co_set;
703     int i;
704    
705 greg 2.1 if (mcat && mcat_last &&
706     !(tmp = rmx_alloc(1, res->imx.ncols, res->rmp->ncomp)))
707     goto memerror;
708     /* figure out what the user set */
709     co_set = fundefined("co");
710     if (!co_set)
711     co_set = -vardefined("co");
712     if (!co_set & (in_ncomp == 3) && vardefined("ro") &&
713     vardefined("go") && vardefined("bo")) {
714     scompile("co(p)=select(p,ro,go,bo)", NULL, 0);
715     co_set = 1;
716     }
717     if (co_set) { /* set if user wants, didn't set */
718     set_r = varlookup("r") != NULL && !vardefined("r");
719     set_c = varlookup("c") != NULL && !vardefined("c");
720     } else /* save a little time */
721     set_r = set_c = 0;
722 greg 2.13
723 greg 2.16 /* read/process row-by-row */
724 greg 2.10 for (cur_row = row0; (in_nrows <= 0) | (cur_row < in_nrows); cur_row += rstep) {
725 greg 2.1 RMATRIX *mres = NULL;
726 greg 2.14 for (i = 0; i < nmats; i++)
727 greg 2.1 if (!rmx_load_row(mop[i].imx.mtx, &mop[i].imx, mop[i].infp)) {
728 greg 2.5 if (cur_row > in_nrows) /* unknown #input rows? */
729 greg 2.10 break;
730 greg 2.1 fprintf(stderr, "%s: read error at row %d\n",
731     mop[i].inspec, cur_row);
732     return(0);
733     }
734 greg 2.10 if (i < nmats)
735     break;
736     for (i = 0; i < nmats; i++)
737 greg 2.1 if (!apply_op(mop[i].rmp, &mop[i].imx, &mop[i].preop))
738     return(0);
739     if (set_r) varset("r", '=', cur_row);
740     for (cur_col = 0; cur_col < in_ncols; cur_col++) {
741     if (set_c) varset("c", '=', cur_col);
742     for (cur_chan = 0; cur_chan < in_ncomp; cur_chan++) {
743     const int ndx = cur_col*in_ncomp + cur_chan;
744     eclock++;
745     if (!co_set) { /* just summing elements? */
746     res->imx.mtx[ndx] = 0;
747     for (i = nmats; i--; )
748     res->imx.mtx[ndx] += mop[i].rmp->mtx[ndx];
749     } else if (co_set > 0) {
750     double dchan = cur_chan+1;
751     res->imx.mtx[ndx] = funvalue("co", 1, &dchan);
752     } else
753     res->imx.mtx[ndx] = varvalue("co");
754     }
755     } /* final conversions */
756     if (!mcat) {
757     if (!apply_op(res->rmp, &res->imx, &res->preop))
758     return(0);
759     } else if (mcat_last) {
760     if (!apply_op(tmp, &res->imx, &res->preop))
761     return(0);
762     mres = rmx_multiply(tmp, mcat);
763     if (!mres)
764     goto multerror;
765     if (!rmx_transfer_data(res->rmp, mres, 0))
766     return(0);
767     } else /* mcat && !mcat_last */ {
768     mres = rmx_multiply(&res->imx, mcat);
769     if (!mres)
770     goto multerror;
771     if (!apply_op(res->rmp, mres, &res->preop))
772     return(0);
773     }
774     rmx_free(mres); mres = NULL;
775 greg 2.13 if (inchild >= 0) { /* children share stdout */
776 greg 2.16 i = getc(stdin); /* signals it's our turn */
777     if (i != EOF) ungetc(i, stdin);
778 greg 2.13 }
779 greg 2.1 if (!rmx_write_data(res->rmp->mtx, res->rmp->ncomp,
780 greg 2.10 res->rmp->ncols, res->rmp->dtype, stdout))
781     return(0);
782 greg 2.16 if (inchild >= 0 && fflush(stdout) == EOF)
783     return(0);
784 greg 2.1 }
785 greg 2.10 return(inchild >= 0 || fflush(stdout) != EOF);
786 greg 2.1 memerror:
787     fputs("Out of buffer space in combine_input()\n", stderr);
788     return(0);
789     multerror:
790     fputs("Unexpected matrix multiply error in combine_input()\n", stderr);
791     return(0);
792     }
793    
794     static int
795     get_factors(double da[], int n, char *av[])
796     {
797     int ac;
798    
799     for (ac = 0; ac < n && isflt(av[ac]); ac++)
800     da[ac] = atof(av[ac]);
801     return(ac);
802     }
803    
804     static void
805     resize_inparr(int n2alloc)
806     {
807     int i;
808    
809 greg 2.6 if (n2alloc == nall)
810     return;
811     for (i = nall; i > n2alloc; i--) {
812 greg 2.1 rmx_reset(&mop[i].imx);
813     if (mop[i].rmp != &mop[i].imx)
814     rmx_free(mop[i].rmp);
815     }
816     mop = (ROPMAT *)realloc(mop, n2alloc*sizeof(ROPMAT));
817     if (mop == NULL) {
818     fputs("Out of memory in resize_inparr()\n", stderr);
819     exit(1);
820     }
821 greg 2.6 if (n2alloc > nall)
822     memset(mop+nall, 0, (n2alloc-nall)*sizeof(ROPMAT));
823 greg 2.1 nall = n2alloc;
824     }
825    
826     /* Load one or more matrices and operate on them, sending results to stdout */
827     int
828     main(int argc, char *argv[])
829     {
830    
831     int outfmt = DTfromHeader;
832     const char *defCsym = NULL;
833     int echoheader = 1;
834     int stdin_used = 0;
835 greg 2.10 int nproc = 1;
836 greg 2.1 const char *mcat_spec = NULL;
837     int n2comp = 0;
838     uby8 comp_ndx[128];
839     int i;
840     /* get starting input array */
841     mop = (ROPMAT *)calloc(nall=2, sizeof(ROPMAT));
842     /* get options and arguments */
843     for (i = 1; i < argc; i++)
844     if (argv[i][0] != '-' || !argv[i][1]) {
845     if (argv[i][0] == '-') {
846     if (stdin_used++) goto stdin_error;
847     mop[nmats].inspec = stdin_name;
848     } else
849     mop[nmats].inspec = argv[i];
850     if (!mop[nmats].preop.csym)
851     mop[nmats].preop.csym = defCsym;
852     if (++nmats >= nall)
853     resize_inparr(nmats + (nmats>>2) + 2);
854     } else {
855     int n = argc-1 - i;
856     switch (argv[i][1]) { /* get option */
857     case 'w':
858     nowarn = !nowarn;
859     break;
860     case 'h':
861     echoheader = !echoheader;
862     break;
863 greg 2.10 case 'n':
864     nproc = atoi(argv[++i]);
865     if (nproc <= 0)
866     goto userr;
867     break;
868 greg 2.1 case 'e':
869     if (!n) goto userr;
870     comp_ndx[n2comp++] = i++;
871     break;
872     case 'f':
873     switch (argv[i][2]) {
874     case '\0':
875     if (!n) goto userr;
876     comp_ndx[n2comp++] = i++;
877     break;
878     case 'd':
879     outfmt = DTdouble;
880     break;
881     case 'f':
882     outfmt = DTfloat;
883     break;
884     case 'a':
885     outfmt = DTascii;
886     break;
887     case 'c':
888     outfmt = DTrgbe;
889     break;
890     default:
891     goto userr;
892     }
893     break;
894     case 's':
895     if (n > MAXCOMP) n = MAXCOMP;
896     i += mop[nmats].preop.nsf =
897     get_factors(mop[nmats].preop.sca,
898     n, argv+i+1);
899     if (mop[nmats].preop.nsf <= 0) {
900     fprintf(stderr, "%s: -s missing arguments\n",
901     argv[0]);
902     goto userr;
903     }
904     break;
905     case 'C':
906     if (!n || isflt(argv[i+1]))
907     goto userr;
908     defCsym = mop[nmats].preop.csym = argv[++i];
909     mop[nmats].preop.clen = 0;
910     mcat_last = 0;
911     break;
912     case 'c':
913     if (n && !isflt(argv[i+1])) {
914     mop[nmats].preop.csym = argv[++i];
915     mop[nmats].preop.clen = 0;
916 greg 2.15 mcat_last = 0;
917 greg 2.1 break;
918     }
919     if (n > MAXCOMP*MAXCOMP) n = MAXCOMP*MAXCOMP;
920     i += mop[nmats].preop.clen =
921     get_factors(mop[nmats].preop.cmat,
922     n, argv+i+1);
923     if (mop[nmats].preop.clen <= 0) {
924     fprintf(stderr, "%s: -c missing arguments\n",
925     argv[0]);
926     goto userr;
927     }
928     mop[nmats].preop.csym = NULL;
929     mcat_last = 0;
930     break;
931     case 'm':
932     if (!n) goto userr;
933     if (argv[++i][0] == '-' && !argv[i][1]) {
934     if (stdin_used++) goto stdin_error;
935     mcat_spec = stdin_name;
936     } else
937     mcat_spec = argv[i];
938     mcat_last = 1;
939     break;
940     default:
941     fprintf(stderr, "%s: unknown option '%s'\n",
942     argv[0], argv[i]);
943     goto userr;
944     }
945     }
946     if (!nmats) {
947     fprintf(stderr, "%s: need at least one input matrix\n", argv[0]);
948     goto userr;
949     }
950     resize_inparr(nmats+1); /* extra matrix at end for result */
951     mop[nmats].inspec = "trailing_ops";
952     /* load final concatenation matrix */
953     if (mcat_spec && !(mcat = rmx_load(mcat_spec, RMPnone))) {
954     fprintf(stderr, "%s: error loading concatenation matrix: %s\n",
955     argv[0], mcat_spec);
956     return(1);
957     }
958     /* get/check inputs, set constants */
959     if (!initialize(&mop[nmats].imx))
960     return(1);
961    
962     for (i = 0; i < n2comp; i++) /* user .cal files and expressions */
963     if (argv[comp_ndx[i]][1] == 'f') {
964     char *fpath = getpath(argv[comp_ndx[i]+1],
965     getrlibpath(), 0);
966     if (fpath == NULL) {
967     fprintf(stderr, "%s: cannot find file '%s'\n",
968     argv[0], argv[comp_ndx[i]+1]);
969     return(1);
970     }
971     fcompile(fpath);
972     } else /* (argv[comp_ndx[i]][1] == 'e') */
973     scompile(argv[comp_ndx[i]+1], NULL, 0);
974    
975     /* get trailing color transform */
976     if (!get_component_xfm(&mop[nmats]))
977     return(1);
978     /* adjust output dimensions and #components */
979     if (mcat) {
980     if (mop[nmats].imx.ncols != mcat->nrows) {
981     fprintf(stderr,
982     "%s: number of input columns does not match number of rows in '%s'\n",
983     argv[0], mcat_spec);
984     return(1);
985     }
986     if (mcat->ncomp != (mcat_last ? mop[nmats].rmp->ncomp : mop[nmats].imx.ncomp)) {
987     fprintf(stderr,
988     "%s: number of components does not match those in '%s'\n",
989     argv[0], mcat_spec);
990     return(1);
991     }
992     if (!split_input(&mop[nmats]))
993     return(1);
994     mop[nmats].rmp->ncols = mcat->ncols;
995     }
996     newheader("RADIANCE", stdout); /* write output header */
997     if (echoheader)
998     output_headinfo(stdout);
999     printargs(argc, argv, stdout);
1000     fputnow(stdout);
1001     mop[nmats].rmp->dtype = rmx_write_header(mop[nmats].rmp, outfmt, stdout);
1002     if (!mop[nmats].rmp->dtype) {
1003     fprintf(stderr, "%s: unsupported output format\n", argv[0]);
1004     return(1);
1005     }
1006 greg 2.7 doptimize(1); /* optimize definitions */
1007 greg 2.11 if (spawned_children(nproc)) /* running in parent process? */
1008 greg 2.10 return(parent_loop() ? 0 : 1);
1009 greg 2.1 /* process & write rows */
1010 greg 2.10 return(combine_input() ? 0 : 1);
1011 greg 2.1 stdin_error:
1012     fprintf(stderr, "%s: %s used for more than one input\n",
1013     argv[0], stdin_name);
1014     return(1);
1015     userr:
1016     fprintf(stderr,
1017 greg 2.10 "Usage: %s [-h][-f{adfc}][-n nproc][-e expr][-f file][-s sf .. | -c ce ..] m1 .. -m mcat > mres\n",
1018 greg 2.1 argv[0]);
1019     return(1);
1020     }