ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/util/rmtxop.c
Revision: 2.32
Committed: Tue Dec 19 16:09:20 2023 UTC (3 months, 1 week ago) by greg
Content type: text/plain
Branch: MAIN
CVS Tags: HEAD
Changes since 2.31: +3 -3 lines
Log Message:
fix(rmtxop,rcomb): Added checks for "xyz" and "rgb" format output

File Contents

# User Rev Content
1 greg 2.1 #ifndef lint
2 greg 2.32 static const char RCSid[] = "$Id: rmtxop.c,v 2.31 2023/12/19 00:39:03 greg Exp $";
3 greg 2.1 #endif
4     /*
5     * General component matrix operations.
6     */
7    
8 greg 2.11 #include <errno.h>
9 greg 2.1 #include "rtio.h"
10     #include "resolu.h"
11     #include "rmatrix.h"
12 greg 2.10 #include "platform.h"
13 greg 2.1
14 greg 2.22 #define MAXCOMP MAXCSAMP /* #components we support */
15 greg 2.1
16 greg 2.22 /* Unary matrix operation(s) */
17 greg 2.1 typedef struct {
18 greg 2.27 double cmat[MAXCOMP*MAXCOMP]; /* component transformation */
19 greg 2.1 double sca[MAXCOMP]; /* scalar coefficients */
20 greg 2.27 const char *csym; /* symbolic coefs or file */
21     short clen; /* number of coefficients */
22 greg 2.13 short nsf; /* number of scalars */
23 greg 2.27 short transpose; /* do transpose? */
24 greg 2.13 } RUNARYOP;
25    
26 greg 2.22 /* Matrix input source and requested operation(s) */
27 greg 2.13 typedef struct {
28     const char *inspec; /* input specification */
29 greg 2.18 RMPref rmp; /* matrix preference */
30 greg 2.13 RUNARYOP preop; /* unary operation(s) */
31     RMATRIX *mtx; /* original matrix if loaded */
32     int binop; /* binary op with next (or 0) */
33     } ROPMAT;
34 greg 2.1
35     int verbose = 0; /* verbose reporting? */
36    
37 greg 2.13 /* Load matrix */
38     static int
39     loadmatrix(ROPMAT *rop)
40 greg 2.1 {
41 greg 2.20 if (rop->mtx != NULL) /* already loaded? */
42 greg 2.13 return(0);
43    
44 greg 2.18 rop->mtx = rmx_load(rop->inspec, rop->rmp);
45 greg 2.23
46     return(!rop->mtx ? -1 : 1);
47 greg 2.1 }
48    
49 greg 2.27 static int checksymbolic(ROPMAT *rop);
50    
51     /* Check/set transform based on a reference input file */
52     static int
53     checkreffile(ROPMAT *rop)
54     {
55     static const char *curRF = NULL;
56     static RMATRIX refm;
57     const int nc = rop->mtx->ncomp;
58     int i;
59    
60     if (!curRF || strcmp(rop->preop.csym, curRF)) {
61     FILE *fp = fopen(rop->preop.csym, "rb");
62     if (!rmx_load_header(&refm, fp)) {
63     fprintf(stderr, "%s: cannot read info header\n",
64     rop->preop.csym);
65     curRF = NULL;
66     if (fp) fclose(fp);
67     return(-1);
68     }
69     fclose(fp);
70     curRF = rop->preop.csym;
71     }
72 greg 2.30 if (refm.ncomp == 3) {
73 greg 2.27 rop->preop.csym = (refm.dtype == DTxyze) ? "XYZ" : "RGB";
74     return(checksymbolic(rop));
75     }
76     if (refm.ncomp == 2) {
77     fprintf(stderr, "%s: cannot convert to 2 components\n",
78     curRF);
79     return(-1);
80     }
81     if (refm.ncomp == 1) {
82     rop->preop.csym = "Y"; /* XXX big assumption */
83     return(checksymbolic(rop));
84     }
85     if (refm.ncomp == nc &&
86     !memcmp(refm.wlpart, rop->mtx->wlpart, sizeof(refm.wlpart)))
87     return(0); /* nothing to do */
88    
89     if ((nc <= 3) | (nc > MAXCSAMP) | (refm.ncomp > MAXCSAMP)) {
90     fprintf(stderr, "%s: cannot resample from %d to %d components\n",
91     curRF, nc, refm.ncomp);
92     return(-1);
93     }
94     rop->preop.clen = refm.ncomp * nc; /* compute spec to ref */
95    
96     for (i = 0; i < nc; i++) {
97     SCOLOR scstim, scresp;
98     int j;
99     memset(scstim, 0, sizeof(COLORV)*nc);
100     scstim[i] = 1.f;
101     convertscolor(scresp, refm.ncomp, refm.wlpart[0], refm.wlpart[3],
102     scstim, nc, rop->mtx->wlpart[0], rop->mtx->wlpart[3]);
103     for (j = refm.ncomp; j-- > 0; )
104     rop->preop.cmat[j*nc + i] = scresp[j];
105     }
106 greg 2.29 memcpy(rop->mtx->wlpart, refm.wlpart, sizeof(rop->mtx->wlpart));
107 greg 2.27 return(0);
108     }
109    
110 greg 2.22 /* Compute conversion row from spectrum to one channel of RGB */
111     static void
112     rgbrow(ROPMAT *rop, int r, int p)
113     {
114     const int nc = rop->mtx->ncomp;
115     const float * wlp = rop->mtx->wlpart;
116     int i;
117    
118     for (i = nc; i--; ) {
119     int nmEnd = wlp[0] + (wlp[3] - wlp[0])*i/nc;
120     int nmStart = wlp[0] + (wlp[3] - wlp[0])*(i+1)/nc;
121     COLOR crgb;
122     spec_rgb(crgb, nmStart, nmEnd);
123     rop->preop.cmat[r*nc+i] = crgb[p];
124     }
125     }
126    
127     /* Compute conversion row from spectrum to one channel of XYZ */
128     static void
129     xyzrow(ROPMAT *rop, int r, int p)
130     {
131     const int nc = rop->mtx->ncomp;
132     const float * wlp = rop->mtx->wlpart;
133     int i;
134    
135     for (i = nc; i--; ) {
136     int nmEnd = wlp[0] + (wlp[3] - wlp[0])*i/nc;
137     int nmStart = wlp[0] + (wlp[3] - wlp[0])*(i+1)/nc;
138     COLOR cxyz;
139     spec_cie(cxyz, nmStart, nmEnd);
140     rop->preop.cmat[r*nc+i] = cxyz[p];
141     }
142     }
143    
144     /* Use the spectral sensitivity function to compute matrix coefficients */
145     static void
146     sensrow(ROPMAT *rop, int r, double (*sf)(SCOLOR sc, int ncs, const float wlpt[4]))
147     {
148     const int nc = rop->mtx->ncomp;
149     int i;
150    
151     for (i = nc; i--; ) {
152     SCOLOR sclr;
153 greg 2.27 memset(sclr, 0, sizeof(COLORV)*nc);
154 greg 2.24 sclr[i] = 1.f;
155 greg 2.22 rop->preop.cmat[r*nc+i] = (*sf)(sclr, nc, rop->mtx->wlpart);
156     }
157     }
158    
159     /* Check/set symbolic transform */
160     static int
161     checksymbolic(ROPMAT *rop)
162     {
163     const int nc = rop->mtx->ncomp;
164     const int dt = rop->mtx->dtype;
165 greg 2.31 double cf = 1;
166 greg 2.22 int i, j;
167 greg 2.27 /* check suffix => reference file */
168     if (strchr(rop->preop.csym, '.') > rop->preop.csym)
169     return(checkreffile(rop));
170 greg 2.22
171     if (nc < 3) {
172     fprintf(stderr, "%s: -c '%s' requires at least 3 components\n",
173     rop->inspec, rop->preop.csym);
174     return(-1);
175     }
176     rop->preop.clen = strlen(rop->preop.csym) * nc;
177     if (rop->preop.clen > MAXCOMP*MAXCOMP) {
178     fprintf(stderr, "%s: -c '%s' results in too many components\n",
179     rop->inspec, rop->preop.csym);
180     return(-1);
181     }
182     for (j = 0; rop->preop.csym[j]; j++) {
183     int comp = 0;
184     switch (rop->preop.csym[j]) {
185     case 'B':
186 greg 2.31 case 'b':
187 greg 2.22 ++comp;
188     /* fall through */
189     case 'G':
190 greg 2.31 case 'g':
191 greg 2.22 ++comp;
192     /* fall through */
193     case 'R':
194 greg 2.31 case 'r':
195     if (rop->preop.csym[j] <= 'Z')
196     cf = 1./WHTEFFICACY;
197 greg 2.22 if (dt == DTxyze) {
198     for (i = 3; i--; )
199 greg 2.31 rop->preop.cmat[j*nc+i] = cf*xyz2rgbmat[comp][i];
200 greg 2.22 } else if (nc == 3)
201     rop->preop.cmat[j*nc+comp] = 1.;
202     else
203     rgbrow(rop, j, comp);
204     break;
205     case 'Z':
206 greg 2.31 case 'z':
207 greg 2.22 ++comp;
208     /* fall through */
209     case 'Y':
210 greg 2.31 case 'y':
211 greg 2.22 ++comp;
212     /* fall through */
213     case 'X':
214 greg 2.31 case 'x':
215     if ((rop->preop.csym[j] <= 'Z') & (dt != DTxyze))
216     cf = WHTEFFICACY;
217 greg 2.22 if (dt == DTxyze) {
218     rop->preop.cmat[j*nc+comp] = 1.;
219     } else if (nc == 3) {
220     for (i = 3; i--; )
221     rop->preop.cmat[j*nc+i] =
222     rgb2xyzmat[comp][i];
223     } else if (comp == CIEY)
224     sensrow(rop, j, scolor2photopic);
225     else
226     xyzrow(rop, j, comp);
227    
228 greg 2.31 for (i = nc*(cf != 1); i--; )
229     rop->preop.cmat[j*nc+i] *= cf;
230 greg 2.22 break;
231 greg 2.24 case 'S': /* scotopic (il)luminance */
232 greg 2.31 cf = WHTSCOTOPIC;
233     /* fall through */
234     case 's':
235 greg 2.22 sensrow(rop, j, scolor2scotopic);
236 greg 2.31 for (i = nc*(cf != 1); i--; )
237     rop->preop.cmat[j*nc+i] *= cf;
238 greg 2.22 break;
239 greg 2.24 case 'M': /* melanopic (il)luminance */
240 greg 2.31 cf = WHTMELANOPIC;
241     /* fall through */
242     case 'm':
243 greg 2.22 sensrow(rop, j, scolor2melanopic);
244 greg 2.31 for (i = nc*(cf != 1); i--; )
245     rop->preop.cmat[j*nc+i] *= cf;
246 greg 2.22 break;
247 greg 2.24 case 'A': /* average component */
248 greg 2.31 case 'a':
249 greg 2.24 for (i = nc; i--; )
250 greg 2.25 rop->preop.cmat[j*nc+i] = 1./(double)nc;
251 greg 2.24 break;
252 greg 2.22 default:
253     fprintf(stderr, "%s: -c '%c' unsupported\n",
254     rop->inspec, rop->preop.csym[j]);
255     return(-1);
256     }
257     }
258     /* return recommended output type */
259 greg 2.32 if (!strcasecmp(rop->preop.csym, "XYZ")) {
260 greg 2.22 if (dt <= DTspec)
261     return(DTxyze);
262 greg 2.32 } else if (!strcasecmp(rop->preop.csym, "RGB")) {
263 greg 2.22 if (dt <= DTspec)
264     return(DTrgbe);
265 greg 2.30 } else if (dt == DTspec)
266 greg 2.22 return(DTfloat); /* probably not actual spectrum */
267     return(0);
268     }
269    
270 greg 2.13 /* Get matrix and perform unary operations */
271 greg 2.1 static RMATRIX *
272 greg 2.13 loadop(ROPMAT *rop)
273 greg 2.1 {
274 greg 2.22 int outtype = 0;
275 greg 2.13 RMATRIX *mres;
276 greg 2.22 int i, j;
277 greg 2.1
278 greg 2.13 if (loadmatrix(rop) < 0) /* make sure we're loaded */
279 greg 2.1 return(NULL);
280 greg 2.13
281 greg 2.27 if (rop->preop.csym && /* symbolic transform? */
282 greg 2.22 (outtype = checksymbolic(rop)) < 0)
283     goto failure;
284     if (rop->preop.clen > 0) { /* apply component transform? */
285     if (rop->preop.clen % rop->mtx->ncomp) {
286     fprintf(stderr, "%s: -c must have N x %d coefficients\n",
287     rop->inspec, rop->mtx->ncomp);
288     goto failure;
289     }
290     if (rop->preop.nsf > 0) { /* scale transform, first */
291     if (rop->preop.nsf == 1) {
292     for (i = rop->preop.clen; i--; )
293     rop->preop.cmat[i] *= rop->preop.sca[0];
294 greg 2.28 } else if (rop->preop.nsf*rop->mtx->ncomp != rop->preop.clen) {
295 greg 2.22 fprintf(stderr, "%s: -s must have one or %d factors\n",
296 greg 2.28 rop->inspec,
297     rop->preop.clen/rop->mtx->ncomp);
298 greg 2.22 goto failure;
299     } else {
300 greg 2.28 for (i = rop->preop.nsf; i--; )
301     for (j = rop->mtx->ncomp; j--; )
302     rop->preop.cmat[i*rop->mtx->ncomp+j]
303     *= rop->preop.sca[i];
304 greg 2.22 }
305     }
306     mres = rmx_transform(rop->mtx, rop->preop.clen/rop->mtx->ncomp,
307     rop->preop.cmat);
308     if (mres == NULL) {
309     fprintf(stderr, "%s: matrix transform failed\n",
310     rop->inspec);
311 greg 2.13 goto failure;
312 greg 2.1 }
313 greg 2.22 if (verbose)
314     fprintf(stderr, "%s: applied %d x %d transform%s\n",
315     rop->inspec, mres->ncomp,
316     rop->mtx->ncomp,
317     rop->preop.nsf ? " (* scalar)" : "");
318 greg 2.26 rop->preop.nsf = 0; /* now folded in */
319 greg 2.22 if ((mres->ncomp > 3) & (mres->dtype <= DTspec))
320     outtype = DTfloat; /* probably not actual spectrum */
321     rmx_free(rop->mtx);
322     rop->mtx = mres;
323     }
324     if (rop->preop.nsf > 0) { /* apply scalar(s)? */
325 greg 2.13 if (rop->preop.nsf == 1) {
326     for (i = rop->mtx->ncomp; --i; )
327     rop->preop.sca[i] = rop->preop.sca[0];
328     } else if (rop->preop.nsf != rop->mtx->ncomp) {
329 greg 2.1 fprintf(stderr, "%s: -s must have one or %d factors\n",
330 greg 2.13 rop->inspec, rop->mtx->ncomp);
331     goto failure;
332 greg 2.1 }
333 greg 2.13 if (!rmx_scale(rop->mtx, rop->preop.sca)) {
334     fputs(rop->inspec, stderr);
335 greg 2.1 fputs(": scalar operation failed\n", stderr);
336 greg 2.13 goto failure;
337 greg 2.1 }
338     if (verbose) {
339 greg 2.13 fputs(rop->inspec, stderr);
340 greg 2.1 fputs(": applied scalar (", stderr);
341 greg 2.13 for (i = 0; i < rop->preop.nsf; i++)
342     fprintf(stderr, " %f", rop->preop.sca[i]);
343 greg 2.1 fputs(" )\n", stderr);
344     }
345     }
346 greg 2.22 if (rop->preop.transpose) { /* transpose matrix? */
347 greg 2.13 mres = rmx_transpose(rop->mtx);
348     if (mres == NULL) {
349     fputs(rop->inspec, stderr);
350 greg 2.12 fputs(": transpose failed\n", stderr);
351 greg 2.13 goto failure;
352 greg 2.12 }
353     if (verbose) {
354 greg 2.13 fputs(rop->inspec, stderr);
355 greg 2.12 fputs(": transposed rows and columns\n", stderr);
356     }
357 greg 2.13 rmx_free(rop->mtx);
358     rop->mtx = mres;
359     }
360     mres = rop->mtx;
361     rop->mtx = NULL;
362 greg 2.22 if (outtype)
363     mres->dtype = outtype;
364 greg 2.13 return(mres);
365     failure:
366     rmx_free(rop->mtx);
367     return(rop->mtx = NULL);
368     }
369    
370     /* Execute binary operation, free matrix arguments and return new result */
371     static RMATRIX *
372     binaryop(const char *inspec, RMATRIX *mleft, int op, RMATRIX *mright)
373     {
374     RMATRIX *mres = NULL;
375     int i;
376    
377     if ((mleft == NULL) | (mright == NULL))
378     return(NULL);
379     switch (op) {
380     case '.': /* concatenate */
381 greg 2.16 if (mleft->ncomp != mright->ncomp) {
382     fputs(inspec, stderr);
383     fputs(": # components do not match\n", stderr);
384     } else if (mleft->ncols != mright->nrows) {
385     fputs(inspec, stderr);
386     fputs(": mismatched dimensions\n",
387     stderr);
388     } else
389     mres = rmx_multiply(mleft, mright);
390 greg 2.13 rmx_free(mleft);
391 greg 2.12 rmx_free(mright);
392 greg 2.1 if (mres == NULL) {
393 greg 2.13 fputs(inspec, stderr);
394 greg 2.16 fputs(": concatenation failed\n", stderr);
395 greg 2.1 return(NULL);
396     }
397     if (verbose) {
398 greg 2.13 fputs(inspec, stderr);
399 greg 2.1 fputs(": concatenated matrix\n", stderr);
400     }
401 greg 2.13 break;
402     case '+':
403     if (!rmx_sum(mleft, mright, NULL)) {
404     fputs(inspec, stderr);
405 greg 2.1 fputs(": matrix sum failed\n", stderr);
406 greg 2.13 rmx_free(mleft);
407 greg 2.1 rmx_free(mright);
408     return(NULL);
409     }
410     if (verbose) {
411 greg 2.13 fputs(inspec, stderr);
412 greg 2.1 fputs(": added in matrix\n", stderr);
413     }
414     rmx_free(mright);
415 greg 2.13 mres = mleft;
416     break;
417     case '*':
418     case '/': {
419     const char * tnam = (op == '/') ?
420 greg 2.11 "division" : "multiplication";
421     errno = 0;
422 greg 2.13 if (!rmx_elemult(mleft, mright, (op == '/'))) {
423 greg 2.11 fprintf(stderr, "%s: element-wise %s failed\n",
424 greg 2.13 inspec, tnam);
425     rmx_free(mleft);
426 greg 2.11 rmx_free(mright);
427     return(NULL);
428     }
429     if (errno)
430     fprintf(stderr,
431     "%s: warning - error during element-wise %s\n",
432 greg 2.13 inspec, tnam);
433 greg 2.11 else if (verbose)
434 greg 2.13 fprintf(stderr, "%s: element-wise %s\n", inspec, tnam);
435 greg 2.11 rmx_free(mright);
436 greg 2.13 mres = mleft;
437     } break;
438     default:
439     fprintf(stderr, "%s: unknown operation '%c'\n", inspec, op);
440     rmx_free(mleft);
441 greg 2.1 rmx_free(mright);
442     return(NULL);
443     }
444 greg 2.13 return(mres);
445     }
446    
447     /* Perform matrix operations from left to right */
448     static RMATRIX *
449     op_left2right(ROPMAT *mop)
450     {
451     RMATRIX *mleft = loadop(mop);
452    
453     while (mop->binop) {
454     if (mleft == NULL)
455     break;
456     mleft = binaryop(mop[1].inspec,
457     mleft, mop->binop, loadop(mop+1));
458     mop++;
459     }
460 greg 2.1 return(mleft);
461     }
462    
463 greg 2.13 /* Perform matrix operations from right to left */
464     static RMATRIX *
465     op_right2left(ROPMAT *mop)
466     {
467     RMATRIX *mright;
468     int rpos = 0;
469     /* find end of list */
470     while (mop[rpos].binop)
471 greg 2.14 if (mop[rpos++].binop != '.') {
472     fputs(
473     "Right-to-left evaluation only for matrix multiplication!\n",
474     stderr);
475     return(NULL);
476     }
477 greg 2.13 mright = loadop(mop+rpos);
478     while (rpos-- > 0) {
479     if (mright == NULL)
480     break;
481 greg 2.20 mright = binaryop(mop[rpos+1].inspec,
482 greg 2.13 loadop(mop+rpos), mop[rpos].binop, mright);
483     }
484     return(mright);
485     }
486    
487     #define t_nrows(mop) ((mop)->preop.transpose ? (mop)->mtx->ncols \
488     : (mop)->mtx->nrows)
489     #define t_ncols(mop) ((mop)->preop.transpose ? (mop)->mtx->nrows \
490     : (mop)->mtx->ncols)
491    
492     /* Should we prefer concatenating from rightmost matrix towards left? */
493     static int
494     prefer_right2left(ROPMAT *mop)
495     {
496     int mri = 0;
497    
498     while (mop[mri].binop) /* find rightmost matrix */
499     if (mop[mri++].binop != '.')
500     return(0); /* pre-empt reversal for other ops */
501    
502     if (mri <= 1)
503     return(0); /* won't matter */
504    
505     if (loadmatrix(mop+mri) < 0) /* load rightmost cat */
506     return(1); /* fail will bail in a moment */
507    
508     if (t_ncols(mop+mri) == 1)
509     return(1); /* definitely better R->L */
510    
511     if (t_ncols(mop+mri) >= t_nrows(mop+mri))
512     return(0); /* ...probably worse */
513    
514     if (loadmatrix(mop) < 0) /* load leftmost */
515     return(0); /* fail will bail in a moment */
516    
517     return(t_ncols(mop+mri) < t_nrows(mop));
518     }
519    
520 greg 2.1 static int
521     get_factors(double da[], int n, char *av[])
522     {
523     int ac;
524    
525     for (ac = 0; ac < n && isflt(av[ac]); ac++)
526     da[ac] = atof(av[ac]);
527     return(ac);
528     }
529    
530 greg 2.13 static ROPMAT *
531 greg 2.26 resize_moparr(ROPMAT *mop, int n2alloc)
532 greg 2.13 {
533     int nmats = 0;
534 greg 2.26 int i;
535 greg 2.13
536     while (mop[nmats++].binop)
537     ;
538 greg 2.26 for (i = nmats; i > n2alloc; i--)
539     rmx_free(mop[i].mtx);
540 greg 2.13 mop = (ROPMAT *)realloc(mop, n2alloc*sizeof(ROPMAT));
541     if (mop == NULL) {
542 greg 2.26 fputs("Out of memory in resize_moparr()\n", stderr);
543 greg 2.13 exit(1);
544     }
545     if (n2alloc > nmats)
546     memset(mop+nmats, 0, (n2alloc-nmats)*sizeof(ROPMAT));
547     return(mop);
548     }
549    
550 greg 2.1 /* Load one or more matrices and operate on them, sending results to stdout */
551     int
552     main(int argc, char *argv[])
553     {
554 greg 2.27 int outfmt = DTfromHeader;
555     const char *defCsym = NULL;
556     int nall = 2;
557     ROPMAT *mop = (ROPMAT *)calloc(nall, sizeof(ROPMAT));
558     int nmats = 0;
559     RMATRIX *mres = NULL;
560     int stdin_used = 0;
561     int i;
562 greg 2.1 /* get options and arguments */
563 greg 2.13 for (i = 1; i < argc; i++) {
564 greg 2.11 if (argv[i][0] && !argv[i][1] &&
565 greg 2.13 strchr(".+*/", argv[i][0]) != NULL) {
566 greg 2.15 if (!nmats || mop[nmats-1].binop) {
567 greg 2.14 fprintf(stderr,
568 greg 2.15 "%s: missing matrix argument before '%c' operation\n",
569 greg 2.14 argv[0], argv[i][0]);
570 greg 2.13 return(1);
571     }
572 greg 2.15 mop[nmats-1].binop = argv[i][0];
573 greg 2.1 } else if (argv[i][0] != '-' || !argv[i][1]) {
574 greg 2.13 if (argv[i][0] == '-') {
575     if (stdin_used++) {
576     fprintf(stderr,
577     "%s: standard input used for more than one matrix\n",
578     argv[0]);
579     return(1);
580     }
581     mop[nmats].inspec = stdin_name;
582     } else
583     mop[nmats].inspec = argv[i];
584 greg 2.27 if (!mop[nmats].preop.csym)
585     mop[nmats].preop.csym = defCsym;
586 greg 2.13 if (nmats > 0 && !mop[nmats-1].binop)
587     mop[nmats-1].binop = '.';
588     nmats++;
589 greg 2.1 } else {
590     int n = argc-1 - i;
591     switch (argv[i][1]) { /* get option */
592     case 'v':
593 greg 2.14 verbose++;
594 greg 2.1 break;
595     case 'f':
596     switch (argv[i][2]) {
597     case 'd':
598     outfmt = DTdouble;
599     break;
600     case 'f':
601     outfmt = DTfloat;
602     break;
603     case 'a':
604     outfmt = DTascii;
605     break;
606     case 'c':
607     outfmt = DTrgbe;
608     break;
609     default:
610     goto userr;
611     }
612     break;
613     case 't':
614 greg 2.13 mop[nmats].preop.transpose = 1;
615 greg 2.1 break;
616     case 's':
617     if (n > MAXCOMP) n = MAXCOMP;
618 greg 2.13 i += mop[nmats].preop.nsf =
619     get_factors(mop[nmats].preop.sca,
620     n, argv+i+1);
621 greg 2.22 if (mop[nmats].preop.nsf <= 0) {
622     fprintf(stderr, "%s: -s missing arguments\n",
623     argv[0]);
624     goto userr;
625     }
626 greg 2.1 break;
627 greg 2.27 case 'C':
628     if (!n || isflt(argv[i+1]))
629     goto userr;
630     defCsym = mop[nmats].preop.csym = argv[++i];
631     mop[nmats].preop.clen = 0;
632     break;
633 greg 2.1 case 'c':
634 greg 2.27 if (n && !isflt(argv[i+1])) {
635     mop[nmats].preop.csym = argv[++i];
636 greg 2.22 mop[nmats].preop.clen = 0;
637     break;
638     }
639 greg 2.1 if (n > MAXCOMP*MAXCOMP) n = MAXCOMP*MAXCOMP;
640 greg 2.13 i += mop[nmats].preop.clen =
641     get_factors(mop[nmats].preop.cmat,
642     n, argv+i+1);
643 greg 2.22 if (mop[nmats].preop.clen <= 0) {
644     fprintf(stderr, "%s: -c missing arguments\n",
645     argv[0]);
646     goto userr;
647     }
648 greg 2.27 mop[nmats].preop.csym = NULL;
649 greg 2.1 break;
650 greg 2.18 case 'r':
651     if (argv[i][2] == 'f')
652     mop[nmats].rmp = RMPreflF;
653     else if (argv[i][2] == 'b')
654     mop[nmats].rmp = RMPreflB;
655     else
656     goto userr;
657     break;
658 greg 2.1 default:
659     fprintf(stderr, "%s: unknown operation '%s'\n",
660     argv[0], argv[i]);
661     goto userr;
662     }
663     }
664 greg 2.14 if (nmats >= nall)
665 greg 2.26 mop = resize_moparr(mop, nall += 2);
666 greg 2.13 }
667     if (mop[0].inspec == NULL) /* nothing to do? */
668 greg 2.1 goto userr;
669 greg 2.15 if (mop[nmats-1].binop) {
670     fprintf(stderr,
671     "%s: missing matrix argument after '%c' operation\n",
672     argv[0], mop[nmats-1].binop);
673     return(1);
674     }
675 greg 2.13 /* favor quicker concatenation */
676 greg 2.14 mop[nmats].mtx = prefer_right2left(mop) ? op_right2left(mop)
677     : op_left2right(mop);
678     if (mop[nmats].mtx == NULL)
679     return(1);
680     /* apply trailing unary operations */
681     mop[nmats].inspec = "trailing_ops";
682     mres = loadop(mop+nmats);
683     if (mres == NULL)
684 greg 2.13 return(1);
685 greg 2.22 if (outfmt == DTfromHeader) /* check data type */
686 greg 2.6 outfmt = mres->dtype;
687 greg 2.22 if (outfmt == DTrgbe) {
688     if (mres->ncomp > 3)
689     outfmt = DTspec;
690     else if (mres->dtype == DTxyze)
691     outfmt = DTxyze;
692     }
693 greg 2.26 newheader("RADIANCE", stdout); /* write result to stdout */
694 greg 2.1 printargs(argc, argv, stdout);
695 greg 2.23 return(rmx_write(mres, outfmt, stdout) ? 0 : 1);
696 greg 2.1 userr:
697     fprintf(stderr,
698 greg 2.26 "Usage: %s [-v][-f{adfc}][-t][-s sf .. | -c ce ..][-rf|-rb] m1 [.+*/] .. > mres\n",
699 greg 2.1 argv[0]);
700     return(1);
701     }