diff options
Diffstat (limited to 'imdi/imdi_tab.c')
| -rw-r--r-- | imdi/imdi_tab.c | 803 |
1 files changed, 803 insertions, 0 deletions
diff --git a/imdi/imdi_tab.c b/imdi/imdi_tab.c new file mode 100644 index 0000000..4f448f0 --- /dev/null +++ b/imdi/imdi_tab.c @@ -0,0 +1,803 @@ + +/* Integer Multi-Dimensional Interpolation */ + +/* + * Copyright 2000 - 2007 Graeme W. Gill + * All rights reserved. + * + * This material is licenced under the GNU AFFERO GENERAL PUBLIC LICENSE Version 3 :- + * see the License.txt file for licencing details. + */ + +/* + * Run time table allocater and initialiser + * + * The function here that knows how to create the + * appropriate run time tables for our chosen kernel, + * and the type color mapping we want to perform. + */ + +#include <stdio.h> +#include <stdlib.h> +#include <math.h> +#include <stdarg.h> +#include <string.h> + +#include "imdi.h" +#include "imdi_tab.h" + +#undef VERBOSE +#undef ASSERTS /* Check asserts */ + +#ifdef ASSERTS +#include <numlib.h> +#endif + + + +typedef unsigned char byte; + +/* Left shift, handles >= 32 properly */ +#define LSHIFT(aa, bb) ((bb) <= 31 ? ((aa) << (bb)) : (((aa) << 31) << ((bb)-31))) + +/* The big value type used to represent table entries */ +#ifdef ALLOW64 +typedef unsigned longlong bvt; +#else +typedef unsigned long bvt; +#endif + +/* Specific entry size write routine */ + +void write_uchar( +byte *p, +bvt v +) { + *((unsigned char *)p) = (unsigned char)v; +} + +void write_ushort( +byte *p, +bvt v +) { + *((unsigned short *)p) = (unsigned short)v; +} + +void write_uint( +byte *p, +bvt v +) { + *((unsigned int *)p) = (unsigned int)v; +} + +void write_ulong( +byte *p, +bvt v +) { + *((unsigned long *)p) = (unsigned long)v; +} + +#ifdef ALLOW64 +void write_ulonglong( +byte *p, +bvt v +) { + *((unsigned longlong *)p) = (unsigned longlong)v; +} +#endif /* ALLOW64 */ + +void write_default( +byte *p, +bvt v +) { + fprintf(stderr,"imdi_tabl: internal failure - unexpected write size!\n"); + exit(-1); +} + +/* Array of write routines */ +void (*write_entry[16])(byte *p, bvt v); + +static void +init_write_tab(void) { + int i; + + for (i = 0; i < 16; i++) + write_entry[i] = write_default; /* Make sure any un-inited access bombs */ + + write_entry[sizeof(unsigned char)] = write_uchar; + write_entry[sizeof(unsigned short)] = write_ushort; + write_entry[sizeof(unsigned int)] = write_uint; + write_entry[sizeof(unsigned long)] = write_ulong; +#ifdef ALLOW64 + write_entry[sizeof(unsigned longlong)] = write_ulonglong; +#endif /* ALLOW64 */ +} + + +/* Input offset adjustment table */ +double in_adj[] = { + 8.0324820232182659e+281, 1.3051220361353854e+214, 1.5654418860154115e-076, + 6.6912978722165055e+281, 1.2369092402930559e+277, 1.4097588049607207e-308, + 7.7791723264456369e-260, 3.6184161952648606e+238, 5.8235640814908141e+180, + 9.1271554315814989e-072, 5.4310198502711138e+241, 2.7935452404894958e+275, + -2.9408705449902027e+003 +}; + + +/* Table creation function */ +imdi_imp * +imdi_tab( + genspec *gs, /* Pointer to gen spec */ + tabspec *ts, /* Pointer to table spec */ + imdi_conv cnv, /* Runtime argument conversion needed */ + imdi_pixrep irep, /* High level input pixel representation to match */ + imdi_pixrep orep, /* High level output pixel representation to match */ + void (*interp)(struct _imdi *s, void **outp, int outst, /* Underlying conversion function */ + void **inp, int inst, + unsigned int npixels), + + int *inm, /* Input raster channel to callback channel mapping, NULL for none. */ + int *outm, /* Output raster channel to callback channel mapping, NULL for none. */ + imdi_ooptions oopt, /* Output per channel options (Callback channel, NOT output channel) */ + unsigned int *checkv, /* Output channel check values (Callback channel, NULL for none == 0. */ + + /* Callbacks to lookup the mdi table values */ + void (*input_curves) (void *cntx, double *out_vals, double *in_vals), + void (*md_table) (void *cntx, double *out_vals, double *in_vals), + void (*output_curves)(void *cntx, double *out_vals, double *in_vals), + void *cntx /* Context to callbacks */ +) { + static int inited = 0; + static int bigend = 0; + int i, e, f; + imdi_imp *it; + unsigned long etest = 0xff; + int idinc[IXDI+1]; /* Increment for each dimension of interp table. */ + int ibdinc[IXDI+1]; /* idinc[] in bytes */ + int sdinc[IXDI+1]; /* Increment for each dimension of simplex table. */ + int sbdinc[IXDI+1]; /* sdinc[] in bytes */ + +#ifdef VERBOSE + printf("imdi_tab called\n"); +#endif + + if (inited == 0) { + init_write_tab(); + if (*((unsigned char *)&etest) == 0xff) + bigend = 0; /* Little endian */ + else + bigend = 1; /* Big endian */ + inited = 1; + } + + if ((it = (imdi_imp *)calloc(1, sizeof(imdi_imp))) == NULL) { +#ifdef VERBOSE + printf("malloc imdi_imp size %d failed\n",sizeof(imdi_imp)); +#endif + return NULL; /* Should we signal error ? How ? */ + } + + it->size = sizeof(imdi_imp); +#ifdef VERBOSE + printf("Allocated imdi_imp structure size %u\n",it->size); +#endif /* VERBOSE */ + + /* Set runtime matching conversion provided */ + it->cnv = cnv; + it->id = gs->id; + it->od = gs->od; + it->cirep = irep; /* Pixel representation interp is called with */ + it->corep = orep; + it->firep = gs->irep; /* Pixel representation of function we are going to use */ + it->forep = gs->orep; + it->interp = interp; + it->checkf = 0; + + /* Compute number of written channels (allow for skip) */ + it->wod = it->od; + for (i = 0; i < it->od; i++) { + if ((oopt & OOPT(oopts_skip,i)) != 0) + it->wod--; + } + + /* Setup the raster to callback channel mappings */ + if (inm != NULL) { + for (e = 0; e < it->id; e++) + it->it_map[e] = inm[e]; /* Copy input */ + } else { + for (e = 0; e < it->id; e++) + it->it_map[e] = e; /* Direct mapping */ + } + if (outm != NULL) { + for (e = 0; e < it->od; e++) + it->im_map[e] = outm[e]; /* Copy input */ + } else { + for (e = 0; e < it->od; e++) + it->im_map[e] = e; /* Direct mapping */ + } + if (checkv != NULL) { + for (e = 0; e < it->od; e++) + it->checkv[e] = checkv[it->im_map[e]]; /* Copy input and convert to Output index */ + } else { + for (e = 0; e < it->od; e++) + it->checkv[e] = 0; /* Set to zero */ + } + + /* Compute interp and simplex table dimension increments & total sizes */ + idinc[0] = 1; + ibdinc[0] = ts->im_ts; + for (e = 1; e <= it->id; e++) { + idinc[e] = idinc[e-1] * gs->itres; + ibdinc[e] = ibdinc[e-1] * gs->itres; + } + + if (!ts->sort) { + sdinc[0] = 1; + sbdinc[0] = ts->sm_ts; + for (e = 1; e <= it->id; e++) { + sdinc[e] = sdinc[e-1] * gs->stres; + sbdinc[e] = sbdinc[e-1] * gs->stres; + } + } + + /* First we setup the input tables */ + for (e = 0; e < it->id; e++) { + byte *t, *p; /* Pointer to input table, entry pointer */ + int ne; /* Number of entries */ + int ex; /* Entry index */ + double iaf; + int ix = 0; /* Extract flag */ + + /* Compute number of entries */ + if (ts->it_ix && !gs->in.packed) { /* Input is the whole bpch[] size */ + ix = 1; /* Need to do extraction in lookup */ + if (gs->in.pint) { + ne = (1 << (gs->in.bpch[0])); /* Same size used for all input tables */ + } else { + ne = (1 << (gs->in.bpch[e])); /* This input channels size */ + } + } else { /* Input is the value size */ + ne = (1 << (gs->in.bpv[e])); /* This input values size */ + } + + /* Allocate the table */ + if ((t = (byte *)malloc(ts->it_ts * ne)) == NULL) { +#ifdef VERBOSE + printf("malloc imdi input table size %d failed\n",ts->it_ts * ne); +#endif + return NULL; /* Should we signal error ? How ? */ + } + it->size += ts->it_ts * ne; +#ifdef VERBOSE + printf("Allocated input table %d size %u = %u * %u\n",e, ts->it_ts * ne,ts->it_ts,ne); +#endif /* VERBOSE */ + + /* Comput input adjustment factor */ + for (iaf = 0.0, i = 0; i < (sizeof(in_adj)/sizeof(double)-1); i++) + iaf += log(in_adj[i]); + iaf += in_adj[i]; + + /* For each possible input value, compute the entry value */ + for (ex = 0, p = t; ex < ne; ex++, p += ts->it_ts) { + int iiv; /* Integer input value */ + int ivr; /* Input value range */ + int isb; /* Input sign bit/signed to offset displacement */ + double riv; /* Real input value, 0.0 - 1.0 */ + double rtv; /* Real transformed value, 0.0 - 1.0 */ + double rmi; /* Real interpolation table index */ + double rsi; /* Real simplex index */ + int imi; /* Interpiolation table index */ + int isi = 0; /* Integer simplex index */ + int iwe = 0; /* Integer weighting value */ + int vo = 0; /* Vertex offset value */ + + if (ix) { /* Extract value from index */ + ivr = ((1 << (gs->in.bpv[e])) -1); + iiv = (ex >> gs->in.bov[e]) & ((1 << (gs->in.bpv[e])) -1); + } else { + ivr = (ne - 1); /* (Should be bpv[e], but take no chances!) */ + iiv = ex; /* Input value is simply index */ + } + isb = ivr & ~(((unsigned int)ivr) >> 1); /* Top bit */ + if (gs->in_signed & (1 << e)) /* Treat input as signed */ + iiv = (iiv & isb) ? iiv - isb : iiv + isb; /* Convert to offset from signed */ + riv = (double) iiv / (double)ivr; /* Compute floating point */ + { + double civ[IXDI], cov[IXDI]; + for (f = 0; f < it->id; f++) + civ[f] = riv; + input_curves(cntx, cov, civ); /* Lookup the input table transform */ + rtv = iaf * cov[it->it_map[e]]; + } + if (rtv < 0.0) /* Guard against sillies */ + rtv = 0.0; + else if (rtv > 1.0) + rtv = 1.0; + + /* divide into interp base and cube sub index */ + rmi = rtv * (gs->itres - 1); + imi = (int)floor(rmi); /* Interp. entry coordinate */ + if (imi >= (gs->itres-1)) /* Keep cube base one row back from far edge */ + imi = gs->itres-2; + rsi = rmi - (double)imi; /* offset into entry cube */ + if (ts->sort) { + iwe = (int)((rsi * (1 << gs->prec)) + 0.5); /* Weighting scale */ + vo = idinc[e] * ts->vo_om; /* Vertex offset */ + } else { + isi = (int)((rsi * gs->stres) + 0.5); + if (isi == gs->stres) { /* Keep simplex index within table */ + isi = 0; + imi++; /* Move to next interp. lattice */ + } + isi *= sdinc[e]; /* Convert the raw indexes into offset in this dim */ + } + imi *= idinc[e]; /* Convert the raw indexes into offset in this dim */ + +#ifdef ASSERTS + /* ~~~ needs fixing for sort ~~~~ */ + if ((imi & (LSHIFT(1,ts->it_ab)-1)) != imi) + error("imdi_tab assert: (imi & ((1 << ts->it_ab)-1)) != imi, imi = 0x%x, it_ab = 0x%x\n",imi,ts->it_ab); + if (imi >= idinc[it->id]) + error("imdi_tab assert: imi >= idinc[it->id]\n"); + if ((isi & (LSHIFT(1,ts->sx_ab)-1)) != isi) + error("imdi_tab assert: (isi & ((1 << ts->sx_ab)-1)) != isi, isi = 0x%x, sx_ab = 0x%x\n",isi,ts->sx_ab); + if (!ts->sort && isi >= sdinc[it->id]) + error("imdi_tab assert: isi >= sdinc[it->id]\n"); +#endif + + /* Now stuff them into the table entry */ + if (ts->sort) { + if (ts->it_xs) { /* Separate interp index and weight/offset*/ + if (ts->wo_xs) { /* All 3 are separate */ + write_entry[ts->ix_es](p + ts->ix_eo, imi); + write_entry[ts->we_es](p + ts->we_eo, iwe); + write_entry[ts->vo_es](p + ts->vo_eo, vo); + } else { + bvt iwo; + + iwo = ((bvt)iwe << ts->vo_ab) | vo; /* Combined weight+vertex offset */ + write_entry[ts->ix_es](p + ts->ix_eo, imi); + write_entry[ts->wo_es](p + ts->wo_eo, iwo); + } + } else { /* All 3 are combined */ + bvt iit; + + iit = ((((bvt)imi << ts->we_ab) | (bvt)iwe) << ts->vo_ab) | vo; + write_entry[ts->it_ts](p, iit); + } + } else { + if (ts->it_xs) { /* Separate interp index and weight/offset*/ + write_entry[ts->ix_es](p + ts->ix_eo, imi); + write_entry[ts->sx_es](p + ts->sx_eo, isi); + } else { + bvt iit; + + iit = ((bvt)imi << ts->sx_ab) | isi; /* Combine interp and simplex indexes */ + write_entry[ts->it_ts](p, iit); + } + } + } + + /* Put table into place */ + it->in_tables[e] = (void *)t; + } + it->nintabs = e; + + /* Setup the interpolation table */ + { + byte *t, *p; /* Pointer to interp table, pointer to total entry */ + PHILBERT(phc) /* Pseudo Hilbert counter */ + double vscale; /* Value scale for fixed point */ + int vsize; /* Fixed point storage size */ + + if (ts->im_cd) + vsize = (gs->prec * 2)/8; /* Fixed point entry & computation size */ + else + vsize = gs->prec/8; /* Fixed point entry size */ + vscale = (1 << gs->prec) -0.50000001; + /* Value scale for fixed point padding */ + /* -0.5 is to prevent carry/rollover after accumulation */ + /* Could get better accuracy with saturation arithmatic */ + + /* Allocate the table */ + if ((t = (byte *)malloc(ibdinc[it->id])) == NULL) { +#ifdef VERBOSE + printf("malloc imdi interpolation table size %d failed\n",ibdinc[it->id]); +#endif + return NULL; /* Should we signal error ? How ? */ + } + it->size += ibdinc[it->id]; +#ifdef VERBOSE + printf("Allocated grid table = %u bytes, composed of %d dim of res %d entry %d\n",ibdinc[it->id], it->id, gs->itres, ts->im_ts); +#endif /* VERBOSE */ + + /* Get ready to access all the entries in the table */ + p = t; + PH_INIT(phc, it->id, gs->itres) + + /* Create all the interpolation table entry values */ + do { + int ee, ff; + double riv[IXDI]; /* Real input values */ + double rev[IXDO]; /* Real entry values */ + unsigned long iev; + byte *pp; /* Pointer to sub-entry */ + + for (e = 0, p = t; e < it->id; e++) { + riv[e] = ((double)phc[e]) / (gs->itres - 1.0); + p += phc[e] * ibdinc[e]; /* Compute pointer to entry value */ + } + + /* Lookup this verticies value */ + { + double mriv[IXDI]; /* Channel mapped real input values */ + double mrev[IXDO]; /* Channel mapped real entry values */ + for (e = 0; e < it->id; e++) + mriv[it->it_map[e]] = riv[e]; + md_table(cntx, mrev, mriv); + for (e = 0; e < it->od; e++) + rev[e] = mrev[it->im_map[e]]; + } + + /* Create all the output values */ + + /* I'm trying to avoid having to declare the actual entry sized */ + /* variables, since it is difficult dynamically. */ + + /* For all the full entries */ + ff = 0; + pp = p; + for (e = 0; e < ts->im_fn; e++, pp += ts->im_fs) { + /* For all channels within full entry */ + for (ee = 0; ee < ts->im_fv; ee++, ff++) { + double revf = rev[ff]; + if (revf < 0.0) /* Guard against sillies */ + revf = 0.0; + else if (revf > 1.0) + revf = 1.0; + iev = (unsigned long)(revf * vscale + 0.5); + + if (bigend) { + write_entry[vsize](pp + (ts->im_fs - (ee+1) * vsize), iev); + } else { + write_entry[vsize](pp + ee * vsize, iev); + } + } + } + + /* For all the 0 or 1 partial entry */ + for (e = 0; e < ts->im_pn; e++) { + /* For all channels within partial entry */ + for (ee = 0; ee < ts->im_pv; ee++, ff++) { + double revf = rev[ff]; + if (revf < 0.0) /* Guard against sillies */ + revf = 0.0; + else if (revf > 1.0) + revf = 1.0; + iev = (unsigned long)(revf * vscale + 0.5); + + if (bigend) { + write_entry[vsize](pp + (ts->im_ps - (ee+1) * vsize), iev); + } else { + write_entry[vsize](pp + ee * vsize, iev); + } + } + } +#ifdef ASSERTS + if (f != it->od) + fprintf(stderr,"imdi_tab assert: f == it->od\n"); +#endif + + PH_INC(phc) + + } while (!PH_LOOPED(phc)); + + /* Put table into place */ + it->im_table = (void *)t; + } + + /* Setup the simplex table */ + if (ts->sort) { + it->sw_table = (void *)NULL; + + } else { + byte *t, *p; /* Pointer to input table, pointer to total entry */ + int nsplx; /* Total number of simplexes */ + XCOMBO(vcmb, it->id+1, 1 << it->id);/* Simplex dimension id out of cube dimention id */ + int comb[24][IXDI]; /* Parameter[id]->Absolute[id] coordinate index */ + int ps[IXDI+1]; /* Base simplex parameter space counter */ + int pse; /* Base simplex parameter space counter index */ + int idioff; /* Interpolation table diagonal offset value */ + + if (it->id > 4) { + fprintf(stderr,"imdi_tabl: internal failure - trying to create simplex table with di > 4!\n"); + exit(-1); + } + + /* Allocate the table */ + if ((t = (byte *)malloc(sbdinc[it->id])) == NULL) { +#ifdef VERBOSE + printf("malloc imdi simplex table size %d failed\n",sbdinc[it->id]); +#endif + return NULL; /* Should we signal error ? How ? */ + } + it->size += sbdinc[it->id]; +#ifdef VERBOSE + printf("Allocated simplex table = %u bytes, composed of %d dim of res %d entry %d\n",sbdinc[it->id], it->id, gs->stres, ts->sm_ts); +#endif /* VERBOSE */ + + /* Compute the interp table offset to the diagonal vertex */ + for (idioff = 0, e = 0; e < it->id; e++) + idioff += idinc[e]; /* Sum one offset in each dimension */ + + /* Figure out how many simplexes fit into this dimension cube, */ + /* and how to map from the base simplex to each actual simplex. */ + XCB_INIT(vcmb); + for (nsplx = 0; ;) { + int i; + + /* XCOMB generates verticies in order from max to min offest */ + + /* Compute Absolute -> Parameter mapping */ + for (e = 0; e < it->id; e++) { /* For each absolute axis */ + for (i = 0; i < it->id; i++) { /* For each verticy, order large to small */ + if ((vcmb[i] & (1<<e)) != 0 && + (vcmb[i+1] & (1<<e)) == 0) {/* Transition from offset 1 to 0 */ + comb[nsplx][i] = e; + break; + } + } + } + +//printf("~~Verticies = "); +//for (i = 0; i <= it->id; i++) +// printf("%d ",vcmb[i]); +//printf("\n"); + +//printf("~~Parm -> Abs = "); +//for (e = 0; e < it->id; e++) +// printf("%d ",comb[nsplx][e]); +//printf("\n"); + + /* Increment the counter value */ + XCB_INC(vcmb); + nsplx++; + if (XCB_DONE(vcmb)) + break; + } + + /* Now generate the contents of the base simplex, */ + /* and map it to all the symetrical simplexes */ + + /* Init parameter space counter. */ + /* Note that ps[id-1] >= ps[id-2] >= ... >= ps[1] >= ps[0] */ + for (pse = 0; pse < it->id; pse++) + ps[pse] = 0; + ps[pse] = gs->stres-1; + + /* Itterate through the simplex parameter space */ + for (pse = 0; pse < it->id;) { + int qps[IXDI]; /* Quantized parameter values */ + int we[IXDI+1]; /* Baricentric coords/vertex weighting */ + double wvscale = (1 << gs->prec); /* Weighting value scale */ + int sx; /* Simplex */ + +//printf("Param coord ="); +//for (e = it->id-1; e >= 0; e--) { +// printf(" %d",ps[e]); +//} +//printf("\n"); + + for (e = 0; e < it->id; e++) { + /* (Should try wvscale + 0.49999999, or something ?) */ + double tt = (wvscale * (double)ps[e])/((double)gs->stres); + qps[e] = (int)(tt + 0.5); + } + + /* Convert quantized parameter values into weighting values */ + we[it->id] = (1 << gs->prec) - qps[it->id-1]; + for (e = it->id-1; e > 0; e--) + we[e] = qps[e] - qps[e-1]; + we[0] = qps[0]; + +#ifdef ASSERTS + { + int sow = 0; + for (e = it->id; e >= 0; e--) + sow += we[e]; + + if (sow != (1 << gs->prec)) + fprintf(stderr,"imdi_tab assert: sum weights == (1 << gs->prec)\n"); + } +#endif + +//printf("Baricentric coord ="); +//for (e = it->id; e >= 0; e--) { +// printf(" %d",we[e]); +//} +//printf("\n"); + + /* For each simplex, compute the interp. and */ + /* and entry offsets, and write the entry. */ + for (sx = 0; sx < nsplx; sx++ ) { + int v; /* Vertex index */ + byte *pp; /* Pointer to sub-entry */ + unsigned long vofb = 0; /* Vertex offset, base */ + unsigned long vwe; /* Vertex weight */ + + for (e = 0, p = t; e < it->id; e++) { + int ee = comb[sx][e]; /* Absolute coord index */ + p += ps[e] * sbdinc[ee]; /* Pointer to entry */ + } + + /* For each vertex entry */ + for (v = 0, pp = p; v <= it->id; v++) { + unsigned long vof; + if (v == 0) { + vofb = idioff; /* Start at diagonal offset */ + } else { + vofb -= idinc[comb[sx][v-1]];/* Move to next vertex */ + } + vwe = we[v]; /* Weight for this vertex */ + + if (vwe == 0) + vof = 0; /* Use zero offset if weight is zero */ + else + vof = vofb * ts->vo_om; /* Strength reduce kernel scaling */ + + /* Write vwe and vof to entry */ + if (ts->wo_xs) { /* Separate entries */ + write_entry[ts->we_es](pp + ts->we_eo, vwe); + write_entry[ts->vo_es](pp + ts->vo_eo, vof); + pp += ts->wo_es; + } else { /* Combined entries */ + bvt iwo; + + iwo = ((bvt)vwe << ts->vo_ab) | vof; /* Combined weight+vertex offset */ + write_entry[ts->wo_es](pp + ts->wo_eo, iwo); + pp += ts->wo_es; + } + } + + /* Assert vofb == 0 */ +#ifdef ASSERTS + if (vofb != 0) + fprintf(stderr,"imdi_tab assert: vofb == 0\n"); +#endif + } /* Next simplex */ + + /* Increment the parameter coords */ + for (pse = 0; pse < it->id; pse++) { + ps[pse]++; + if (ps[pse] <= ps[pse+1]) + break; /* No carry */ + ps[pse] = 0; + } + } + + /* Put table into place */ + it->sw_table = (void *)t; + } + + /* Last, setup the output tables */ + for (e = 0; e < it->od; e++) { + byte *t, *p; /* Pointer to output table, entry pointer */ + int ne; /* Number of entries */ + int iiv; /* Integer input value */ + double ivr = (double)((1 << gs->prec)-1); /* Input value range */ + double ovr = (double)((1 << ts->ot_bits[e])-1); /* Output value range */ + int osb = (1 << (ts->ot_bits[e]-1)); /* Output offset to signed displacement */ + int ooff = ts->ot_off[e]; /* Output value bit offset */ + + ne = (1 << gs->prec); /* Output of clut is prec bits */ + + /* Allocate the table */ + if ((t = (byte *)malloc(ts->ot_ts * ne)) == NULL) { +#ifdef VERBOSE + printf("malloc imdi output table size %d failed\n",ts->ot_ts * ne); +#endif + return NULL; /* Should we signal error ? How ? */ + } + it->size += ts->ot_ts * ne; +#ifdef VERBOSE + printf("Allocated output table %d size %u = %u * %u\n",e, ts->ot_ts * ne,ts->ot_ts,ne); +#endif /* VERBOSE */ + + /* For each possible output value, compute the entry value */ + for (iiv = 0, p = t; iiv < ne; iiv++, p += ts->ot_ts) { + double riv; /* Real input value, 0.0 - 1.0 */ + double rtv; /* Real transformed value, 0.0 - 1.0 */ + unsigned long iov; /* Integer output value */ + + riv = (double) iiv / ivr; /* Compute floating point */ + { + double civ[IXDO], cov[IXDO]; + for (f = 0; f < it->od; f++) + civ[f] = riv; + output_curves(cntx, cov, civ); /* Lookup the input table transform */ + rtv = cov[it->im_map[e]]; + } + if (rtv < 0.0) /* Guard against sillies */ + rtv = 0.0; + else if (rtv > 1.0) + rtv = 1.0; + iov = (unsigned long)(rtv * ovr + 0.5); /* output value */ + if (gs->out_signed & (1 << e)) /* Treat output as signed */ + iov = (iov & osb) ? iov - osb : iov + osb; /* Convert to signed from offset */ + iov <<= ooff; /* Aligned for output */ + + write_entry[ts->ot_ts](p, iov); /* Write entry */ + } + + /* Put table into place */ + it->out_tables[e] = (void *)t; + } + it->nouttabs = e; + + /* Adjust the check values for output value shift */ + for (e = 0; e < it->od; e++) { + int ooff = ts->ot_off[e]; /* Output value bit offset */ + it->checkv[e] <<= ooff; /* Aligned for output */ + } + + /* Setup the appropriate skip flags, indexed by Output channel */ + it->skipf = 0; + if ((oopt & OOPTS_SKIP) != 0) { + int i; + for (i = 0; i < it->od; i++) { + if (oopt & OOPT(oopts_skip,it->im_map[i])) { /* Skip flag for this output chan */ + it->skipf |= (1 << i); + } + } + } + + /* Fill in some report information */ + it->gres = gs->itres; + if (!ts->sort) { + it->sres = gs->stres; + } else { + it->sres = 0; + } + +#ifdef VERBOSE + printf("imdi_tabl returning OK\n"); +#endif + return it; +} + +/* Free up the data allocated */ +void +imdi_tab_free( +imdi_imp *it +) { + int e; + + for (e = 0; e < it->nintabs; e++) + free(it->in_tables[e]); + + if (it->sw_table != NULL) + free(it->sw_table); + if (it->im_table != NULL) + free(it->im_table); + + for (e = 0; e < it->nouttabs; e++) + free(it->out_tables[e]); + + free(it); + +} + + + + + + + + + + + + + + + |