diff options
Diffstat (limited to 'pkg/utilities/nttools/imtab')
| -rw-r--r-- | pkg/utilities/nttools/imtab/imtab.h | 4 | ||||
| -rw-r--r-- | pkg/utilities/nttools/imtab/imtab.x | 476 | ||||
| -rw-r--r-- | pkg/utilities/nttools/imtab/itbwcs.x | 129 | ||||
| -rw-r--r-- | pkg/utilities/nttools/imtab/mkpkg | 13 | ||||
| -rw-r--r-- | pkg/utilities/nttools/imtab/tabim.x | 176 |
5 files changed, 798 insertions, 0 deletions
diff --git a/pkg/utilities/nttools/imtab/imtab.h b/pkg/utilities/nttools/imtab/imtab.h new file mode 100644 index 00000000..2f280d60 --- /dev/null +++ b/pkg/utilities/nttools/imtab/imtab.h @@ -0,0 +1,4 @@ +define IMTAB_NO_WCS 0 # pixel coordinates will not be written out +define IMTAB_LOGICAL 1 # wcs = logical +define IMTAB_PHYSICAL 2 # wcs = physical +define IMTAB_WORLD 3 # wcs = world diff --git a/pkg/utilities/nttools/imtab/imtab.x b/pkg/utilities/nttools/imtab/imtab.x new file mode 100644 index 00000000..8e8cb324 --- /dev/null +++ b/pkg/utilities/nttools/imtab/imtab.x @@ -0,0 +1,476 @@ +include <error.h> +include <fset.h> # to check whether output is redirected +include <imhdr.h> +include <mach.h> # for MAX_SHORT +include <tbset.h> +include "imtab.h" + +define NCOLS (1 + IM_MAXDIM) # max number of columns to write to the table + +# imtab -- create a table from an image +# This task copies data from an image to a table. Pixel values are +# read from the image line by line and written to a column in increasing +# row number. +# If the table already exists then columns will be added to it; note that +# the column names must not conflict with existing names. +# +# Phil Hodge, 10-Jan-1991 Task created. +# Phil Hodge, 17-Mar-1992 Include text as a valid table type; call pargstr +# to pass cname in case column name already exists. +# Phil Hodge, 16-Apr-1993 Include short datatype. +# Phil Hodge, 28-Sep-1993 Include wcs option for pixel coordinates. +# Phil Hodge, 13-Dec-1993 Slight changes to itb_init because of optimizer +# problems with SGI Fortran. +# Phil Hodge, 8-Jun-1999 Set output to STDOUT if redirected. +# Phil Hodge, 30-Mar-2000 Allow lists of names for input and output. + +procedure imtab() + +pointer input # name of an input image +pointer outlist # names of output tables +pointer outtable # name of an output table +char cname[SZ_COLNAME,NCOLS] # column names +char c_root[SZ_COLNAME] # root for column names for position +pointer wcs # wcs name if c_root != "" +pointer formats # list of formats for pixel coords +pointer ttype # type of output table (if new) +#-- +pointer sp +pointer im # pointer to image descriptors +pointer xps, xpr, xpd # pointer to input data from image +pointer tp # pointer to descriptor for output table +pointer cp[NCOLS] # column descriptors +pointer mw, ct # mwcs pointers + +pointer imt, tnt # pointers for filename templates +int nin, nout # numbers of names in lists +int junk + +long v[IM_MAXDIM] # for call to imgnld +int lcoords[IM_MAXDIM] # "logical" coordinates, copied from v +real ipcoords[IM_MAXDIM] # "logical" coordinates, copied from v +real opcoords[IM_MAXDIM] # "physical" coordinates from ipcoords +double iwcoords[IM_MAXDIM] # "logical" coordinates, copied from v +double owcoords[IM_MAXDIM] # "world" coordinates from iwcoords + +int wcs_type # wcs name as an int +int ax[IM_MAXDIM] # ax[i] is physical axis for logical axis i +real inr[IM_MAXDIM] # copy of input coordinates used by itb_ctranr +double ind[IM_MAXDIM] # copy of input coordinates used by itb_ctrand +int wcsdim # dimension of physical image coord system +int impixtype # data type of image + +int ncols # number of columns to create +int dtype[NCOLS] # data type of each table column +int frow, lrow # row number limits for tbcptd +int row # loop index +int i, j, k # loop indexes +bool done # loop-termination flag +int clgwrd() +int fstati() +int imgnls(), imgnlr(), imgnld() +pointer imtopenp(), tbnopen() +int imtlen(), imtgetim(), tbnlen(), tbnget() +bool streq() + +begin + call smark (sp) + call salloc (input, SZ_FNAME, TY_CHAR) + call salloc (outlist, SZ_LINE, TY_CHAR) + call salloc (outtable, SZ_FNAME, TY_CHAR) + call salloc (wcs, SZ_FNAME, TY_CHAR) + call salloc (formats, SZ_FNAME, TY_CHAR) + call salloc (ttype, SZ_FNAME, TY_CHAR) + + imt = imtopenp ("input") + nin = imtlen (imt) + + if (fstati (STDOUT, F_REDIR) == YES) + call strcpy ("STDOUT", Memc[outlist], SZ_LINE) + else + call clgstr ("outtable", Memc[outlist], SZ_LINE) + tnt = tbnopen (Memc[outlist]) + nout = tbnlen (tnt) + + # Compare the numbers of input and output names. + call itb_names (nin, nout, Memc[outlist]) + + # Get the column names. + call clgstr ("colname", cname[1,1], SZ_COLNAME) # name for data + call clgstr ("pname", c_root, SZ_COLNAME) # root name for position + call xt_stripwhite (c_root) + if (c_root[1] != EOS) { + wcs_type = clgwrd ("wcs", Memc[wcs], SZ_FNAME, + "|logical|physical|world") + call clgstr ("formats", Memc[formats], SZ_FNAME) + } else { + wcs_type = IMTAB_NO_WCS + Memc[formats] = EOS + } + + # What table type should be created? + if (streq (Memc[outlist], "STDOUT")) + call strcpy ("default", Memc[ttype], SZ_FNAME) + else + call clgstr ("tbltype", Memc[ttype], SZ_FNAME) + + # Loop over the list of input images. + while (imtgetim (imt, Memc[input], SZ_FNAME) != EOF) { + + if (nout == 1) + call tbnrew (tnt) + junk = tbnget (tnt, Memc[outtable], SZ_FNAME) + + # Open the input image and the wcs and get the column data types. + call itb_init (Memc[input], wcs_type, im, dtype) + impixtype = IM_PIXTYPE(im) + + # Initialize the wcs. + call itb_wcs_init (im, wcs_type, mw, ct, ax, wcsdim) + call amovkr (1., inr, IM_MAXDIM) + call amovkd (1.d0, ind, IM_MAXDIM) + + # Initialize for reading the image. + do k = 1, IM_MAXDIM + v[k] = 1 + + # Initialize for writing to the table. + row = 1 + frow = 1 + lrow = IM_LEN(im,1) + if (wcs_type == IMTAB_NO_WCS) + ncols = 1 + else if (wcs_type == IMTAB_LOGICAL) + ncols = 1 + IM_NDIM(im) + else if (wcs_type == IMTAB_PHYSICAL || wcs_type == IMTAB_WORLD) + ncols = 1 + wcsdim + + # Open or create the output table. + call itb_table (im, Memc[outtable], wcs_type, Memc[ttype], ncols, + cname, c_root, Memc[formats], dtype, tp, cp) + + # Copy each line of the image into a column of the table. + done = false + while (!done) { + + # Assign pixel index for all but the first axis. + if (wcs_type == IMTAB_LOGICAL) + do j = 2, ncols-1 + lcoords[j] = v[j] + else if (wcs_type == IMTAB_PHYSICAL) + do j = 2, ncols-1 + ipcoords[j] = real (v[j]) + else if (wcs_type == IMTAB_WORLD) + do j = 2, ncols-1 + iwcoords[j] = double (v[j]) + + if (impixtype == TY_SHORT || impixtype == TY_UBYTE) + done = (imgnls (im, xps, v) == EOF) + else if (impixtype == TY_REAL || impixtype == TY_USHORT) + done = (imgnlr (im, xpr, v) == EOF) + else + done = (imgnld (im, xpd, v) == EOF) + + if (!done) { + + # Write the pixel coordinates. + if (wcs_type == IMTAB_LOGICAL) { + do i = 1, IM_LEN(im,1) { # simply write pixel numbers + lcoords[1] = i + call tbrpti (tp, cp[2], lcoords, ncols-1, row) + row = row + 1 + } + } else if (wcs_type == IMTAB_PHYSICAL) { + do i = 1, IM_LEN(im,1) { + ipcoords[1] = real (i) + call itb_ctranr (im, ct, ax, inr, + ipcoords, opcoords, wcsdim) + call tbrptr (tp, cp[2], opcoords, ncols-1, row) + row = row + 1 + } + } else if (wcs_type == IMTAB_WORLD) { + do i = 1, IM_LEN(im,1) { + iwcoords[1] = double (i) + call itb_ctrand (im, ct, ax, ind, + iwcoords, owcoords, wcsdim) + call tbrptd (tp, cp[2], owcoords, ncols-1, row) + row = row + 1 + } + } + + # Copy image line into a portion of a column of the table. + if (impixtype == TY_SHORT || impixtype == TY_UBYTE) + call tbcpts (tp, cp[1], Mems[xps], frow, lrow) + else if (impixtype == TY_REAL || impixtype == TY_USHORT) + call tbcptr (tp, cp[1], Memr[xpr], frow, lrow) + else + call tbcptd (tp, cp[1], Memd[xpd], frow, lrow) + frow = frow + IM_LEN(im,1) + lrow = lrow + IM_LEN(im,1) + } + } + if (mw != NULL) + call mw_close (mw) # close mwcs + call imunmap (im) # close image + call tbtclo (tp) # close table + } + + call imtclose (imt) + call tbnclose (tnt) + + call sfree (sp) +end + +# This routine checks the number of input and output file names. +# The number of names in the input and output lists must be the same, +# unless all the input will be written to the standard output. + +procedure itb_names (nin, nout, outlist) + +int nin # i: number of input image names +int nout # i: number of output table names +char outlist[ARB] # i: output names (to be compared with "STDOUT") +#-- +bool strne() + +begin + if (nin == 0) + call error (1, "no input image specified") + + if (nout == 0) + call error (1, "no output table specified") + + if (nin != nout && strne (outlist, "STDOUT")) { + + if (nin == 1) { + call eprintf ("There is one input image") + } else { + call eprintf ("There are %d input images") + call pargi (nin) + } + if (nout == 1) { + call eprintf (" and one output table;\n") + } else { + call eprintf (" and %d output tables;\n") + call pargi (nout) + } + call error (1, "the lists must have the same length") + } +end + +# itb_init -- get data types and column info +# This routine opens the input image and gets the data type for each column. + +procedure itb_init (input, wcs_type, im, dtype) + +char input[ARB] # i: name of image +int wcs_type # i: type of wcs for pixel coordinates +pointer im # o: imhdr pointer +int dtype[NCOLS] # o: data type of each table column +#-- +int i # loop index +int fill_extra # dummy +pointer immap() + +begin + # Open input image. + im = immap (input, READ_ONLY, NULL) + + # Fewer data types are allowed for tables than for images. + switch (IM_PIXTYPE(im)) { + case TY_UBYTE, TY_SHORT: + dtype[1] = TY_SHORT + case TY_USHORT, TY_INT, TY_LONG: + dtype[1] = TY_INT + case TY_REAL: + dtype[1] = TY_REAL + case TY_DOUBLE: + dtype[1] = TY_DOUBLE + default: + call error (1, "image data type not supported for tables") + } + + # Set the data types of columns for pixel coordinates. + fill_extra = IM_NDIM(im) + 2 + if (wcs_type == IMTAB_NO_WCS) { + do i = 2, NCOLS + dtype[i] = TY_SHORT # ignored + + } else if (wcs_type == IMTAB_LOGICAL) { + # Check the image size to see if we can use TY_SHORT. + do i = 1, IM_NDIM(im) { + if (IM_LEN(im,i) > MAX_SHORT) + dtype[i+1] = TY_INT + else + dtype[i+1] = TY_SHORT + } + do i = fill_extra, NCOLS + dtype[i] = TY_SHORT # ignored + + } else if (wcs_type == IMTAB_PHYSICAL) { + do i = 2, NCOLS + dtype[i] = TY_REAL + + } else if (wcs_type == IMTAB_WORLD) { + do i = 2, NCOLS + dtype[i] = TY_DOUBLE + } +end + +# itb_table -- initialization for output table +# This routine opens the output table (or creates it if it doesn't already +# exist) and creates the columns for the data and the pixel coordinates. + +procedure itb_table (im, outtable, wcs_type, ttype, ncols, + cname, c_root, formats, dtype, tp, cp) + +pointer im # i: imhdr pointer for input image +char outtable[ARB] # i: name of output table +int wcs_type # i: type of wcs for pixel coordinates +char ttype[ARB] # i: table type (e.g. "row") +int ncols # i: total number of columns to write +char cname[SZ_COLNAME,NCOLS] # io: column names +char c_root[ARB] # i: root for column name for pixels +char formats[ARB] # i: user-specified formats for pixels +int dtype[NCOLS] # i: data types of table columns +pointer tp # o: pointer to table descriptor +pointer cp[NCOLS] # o: column descriptors +#-- +char colunits[SZ_COLUNITS,NCOLS] # column units +char colfmt[SZ_COLFMT,NCOLS] # column format +char history[SZ_FNAME] # for history records +int lendat[NCOLS] # one +int nrows +int i +bool new_table # true if the table does not already exist +bool column_conflict # true if column already exists +int ip, ctowrd() +pointer tbtopn() +int tbtacc() + +begin + colunits[1,1] = EOS + colfmt[1,1] = EOS + + # Assign column names. + do i = 2, NCOLS { + call sprintf (cname[1,i], SZ_COLNAME, "%s%d") + call pargstr (c_root) + call pargi (i-1) + } + + # Replace commas with blanks in the user-specified format string. + ip = 1 + while (formats[ip] != EOS) { + if (formats[ip] == ',') + formats[ip] = ' ' + ip = ip + 1 + } + + # Assign print format and units. + ip = 1 + do i = 2, NCOLS { + + if (wcs_type == IMTAB_LOGICAL) { + + # If the user specified a format, use it; otherwise, + # assign a default. + if (ctowrd (formats, ip, colfmt[1,i], SZ_COLFMT) < 1) { + if (dtype[i] == TY_INT) + call strcpy ("%11d", colfmt[1,i], SZ_COLFMT) + else if (dtype[i] == TY_SHORT) + call strcpy ("%5d", colfmt[1,i], SZ_COLFMT) + } + call strcpy ("pixels", colunits[1,i], SZ_COLUNITS) + + } else if (wcs_type == IMTAB_PHYSICAL) { + + if (ctowrd (formats, ip, colfmt[1,i], SZ_COLFMT) < 1) + call strcpy ("%9.3f", colfmt[1,i], SZ_COLFMT) + call strcpy ("pixels", colunits[1,i], SZ_COLUNITS) + + } else if (wcs_type == IMTAB_WORLD) { + + if (ctowrd (formats, ip, colfmt[1,i], SZ_COLFMT) < 1) + colfmt[1,i] = EOS # take the default + colunits[1,i] = EOS # we don't know the units + + } else { + colfmt[1,i] = EOS + colunits[1,i] = EOS + } + } + + do i = 1, NCOLS { + lendat[i] = 1 + cp[i] = NULL + } + + nrows = 1 + do i = 1, IM_NDIM(im) + nrows = nrows * IM_LEN(im,i) + + # Does the table already exist? + new_table = (tbtacc (outtable) == NO) + + if (new_table) { + tp = tbtopn (outtable, NEW_FILE, NULL) + + if (ttype[1] == 'r') { + call tbpset (tp, TBL_WHTYPE, TBL_TYPE_S_ROW) + } else if (ttype[1] == 'c') { + call tbpset (tp, TBL_WHTYPE, TBL_TYPE_S_COL) + call tbpset (tp, TBL_ALLROWS, nrows) + } else if (ttype[1] == 't') { + call tbpset (tp, TBL_WHTYPE, TBL_TYPE_TEXT) + } + } else { + tp = tbtopn (outtable, READ_WRITE, NULL) + } + + # Make sure the columns don't already exist. + column_conflict = false + if ( ! new_table ) { + call tbcfnd (tp, cname, cp, ncols) + do i = 1, ncols { + if (cp[i] != NULL) { + call eprintf ("Column %s already exists.\n") + call pargstr (cname[1,i]) + column_conflict = true + } + } + if (column_conflict) { + call imunmap (im) + call tbtclo (tp) + call error (1, + "new columns in existing table must be unique") + } + } + + # Define the columns. + call tbcdef (tp, cp, cname, colunits, colfmt, dtype, lendat, ncols) + + if (new_table) + call tbtcre (tp) # open the file + + # Write history info. + call strcpy ("Column ", history, SZ_FNAME) + call strcat (cname[1,1], history, SZ_FNAME) # column name for data + call strcat (" from ", history, SZ_FNAME) + call strcat (IM_HDRFILE(im), history, SZ_FNAME) # name of input image + call tbhadt (tp, "history", history) + if (ncols > 1) { + call strcpy ("Column ", history, SZ_FNAME) + call strcat (cname[1,1], history, SZ_FNAME) + if (ncols > 2) + call strcat (", pixel columns ", history, SZ_FNAME) + else + call strcat (", pixel column ", history, SZ_FNAME) + do i = 2, ncols-1 { + call strcat (cname[1,i], history, SZ_FNAME) + call strcat (", ", history, SZ_FNAME) + } + call strcat (cname[1,ncols], history, SZ_FNAME) + call tbhadt (tp, "history", history) + } +end diff --git a/pkg/utilities/nttools/imtab/itbwcs.x b/pkg/utilities/nttools/imtab/itbwcs.x new file mode 100644 index 00000000..4ecea420 --- /dev/null +++ b/pkg/utilities/nttools/imtab/itbwcs.x @@ -0,0 +1,129 @@ +include <imhdr.h> +include <mwset.h> +include "imtab.h" + +# This file contains three routines: itb_wcs_init and the single and +# double precision routines itb_ctranr & itb_ctrand. +# +# Phil Hodge, 30-Sep-1993 Subroutines created. + +# itb_wcs_init -- open wcs, etc +# This routine gets the wcs, turns axis mapping off, and initializes +# the transformation for physical or world coordinates. The dimension +# of the original image and the mapping from logical to physical axis +# numbers are returned for use by itb_ctrand and itb_ctranr. + +procedure itb_wcs_init (im, wcs_type, mw, ct, ax, wcsdim) + +pointer im # i: imhdr pointer for image +int wcs_type # i: wcs type +pointer mw, ct # o: mwcs pointers +int ax[IM_MAXDIM] # o: ax[i] is physical axis for logical axis i +int wcsdim # o: dimension of physical image coord system +#-- +int axno[IM_MAXDIM] # axno[j] is logical axis for physical axis j +int axval[IM_MAXDIM] # axval[j] is value if axno[j] is zero +int ndim # number of "logical" axes +int i, j +pointer mw_openim(), mw_sctran() +int mw_stati() + +begin + if (wcs_type == IMTAB_NO_WCS || wcs_type == IMTAB_LOGICAL) { + mw = NULL + ct = NULL + wcsdim = IM_NDIM(im) + return + } + + # Get the wcs. + mw = mw_openim (im) + + # Set up the transformation. + call mw_seti (mw, MW_USEAXMAP, NO) + if (wcs_type == IMTAB_PHYSICAL) + ct = mw_sctran (mw, "logical", "physical", 0) + else if (wcs_type == IMTAB_WORLD) + ct = mw_sctran (mw, "logical", "world", 0) + wcsdim = mw_stati (mw, MW_NPHYSDIM) + ndim = IM_NDIM(im) + + # Get the logical axis number corresponding to each physical axis. + call mw_gaxmap (mw, axno, axval, wcsdim) + + # Invert axno: get the physical axis number for each logical axis. + do i = 1, ndim # initialize + ax[i] = 0 + do j = 1, wcsdim { + do i = 1, ndim { + if (axno[j] == i) { + ax[i] = j + break + } + } + } + + # Check to be sure each axis was found. + do i = 1, ndim { + if (ax[i] < 1) + call error (1, "itb_mwcs_init: an axis was not found") + } +end + +# itb_ctran -- translate coordinates with axis mapping = NO +# This routine translates "logical" coordinates to "physical" or "world". +# Axis mapping must have been turned off, and the mapping from logical +# to physical axes is given by the array AX: if I is a logical axis +# number, AX[I] is the corresponding physical axis number. Each element +# of the array IN must have been initialized to one by the calling routine. +# Separate single and double precision versions are included. + +procedure itb_ctrand (im, ct, ax, in, incoords, outcoords, wcsdim) + +pointer im # i: imhdr pointer +pointer ct # i: coordinate transformation pointer +int ax[wcsdim] # i: "logical" to "physical" mapping +double in[IM_MAXDIM] # io: copy of incoords but includes axis mapping +double incoords[wcsdim] # i: input "logical" coordinates +double outcoords[wcsdim] # o: output coordinates +int wcsdim # i: length of incoords & outcoords arrays +#-- +int i + +begin + if (ct == NULL) { + call amovd (incoords, outcoords, wcsdim) + return + } + + # Take account of axis mapping; i is the logical axis number. + do i = 1, IM_NDIM(im) + in[ax[i]] = incoords[i] + + call mw_ctrand (ct, in, outcoords, wcsdim) +end + +procedure itb_ctranr (im, ct, ax, in, incoords, outcoords, wcsdim) + +pointer im # i: imhdr pointer +pointer ct # i: coordinate transformation pointer +int ax[wcsdim] # i: "logical" to "physical" mapping +real in[IM_MAXDIM] # io: copy of incoords but includes axis mapping +real incoords[wcsdim] # i: input "logical" coordinates +real outcoords[wcsdim] # o: output coordinates +int wcsdim # i: length of incoords & outcoords arrays +#-- +int i + +begin + if (ct == NULL) { + call amovr (incoords, outcoords, wcsdim) + return + } + + # Take account of axis mapping; i is the logical axis number. + do i = 1, IM_NDIM(im) + in[ax[i]] = incoords[i] + + call mw_ctranr (ct, in, outcoords, wcsdim) +end diff --git a/pkg/utilities/nttools/imtab/mkpkg b/pkg/utilities/nttools/imtab/mkpkg new file mode 100644 index 00000000..6e6431a0 --- /dev/null +++ b/pkg/utilities/nttools/imtab/mkpkg @@ -0,0 +1,13 @@ +# Update the imtab & tabim application code in the ttools package library. +# Author: HODGE, 31-DEC-1990 + +$checkout libpkg.a ../ +$update libpkg.a +$checkin libpkg.a ../ +$exit + +libpkg.a: + imtab.x <error.h> <imhdr.h> <mach.h> <tbset.h> "imtab.h" + itbwcs.x <imhdr.h> <mwset.h> "imtab.h" + tabim.x <imhdr.h> <tbset.h> + ; diff --git a/pkg/utilities/nttools/imtab/tabim.x b/pkg/utilities/nttools/imtab/tabim.x new file mode 100644 index 00000000..8a9e7211 --- /dev/null +++ b/pkg/utilities/nttools/imtab/tabim.x @@ -0,0 +1,176 @@ +include <imhdr.h> +include <fset.h> # to check whether input is redirected +include <tbset.h> + +# tabim -- create an image from one column of a table +# This task copies a column of a table into an image. If the image already +# exists, it will be overwritten; otherwise, a new image will be created. +# For a new image, if the 'ndim' parameter is greater than zero, the size of +# the image will be taken from the parameters 'n1', 'n2', etc. It is the +# user's responsibility to ensure that the product of these values equals +# the number of rows in the table. +# +# Phil Hodge, 12-Oct-1989 Task created +# Phil Hodge, 11-Jan-1991 Allow multi-dimensional output. +# Phil Hodge, 15-May-1998 Check null flag, and replace INDEF with -999. +# Phil Hodge, 8-Jun-1999 Set input to STDIN if redirected. +# Phil Hodge, 30-Mar-2000 Allow lists of names for input and output. + +procedure tabim() + +pointer inlist, outlist # for input and output names +char intable[SZ_FNAME] # name of an input table +char output[SZ_FNAME] # name of an output image +char colname[SZ_COLNAME] # column name +int ndim # dimension of output image +int axlen[IM_MAXDIM] # length of each axis of output image +#-- +pointer sp # stack pointer for scratch space +pointer tp # pointer to descriptor for input table +pointer cp # column descriptor +pointer im # pointer to image descriptor +pointer xp # pointer to output data for image +pointer temp # scratch for parameter name +pointer nullflag # scratch for null flags (ignored) +long v[IM_MAXDIM] # for call to impnld() +int dtype # data type +int npix # number of pixels, accumulated one axis at a time +int nrows # number of rows in table +int nlines # number of lines in image +int frow, lrow # row number limits for tbcgtd +int i, k +int junk +bool new_image # true if the image does not already exist +pointer immap(), tbtopn() +int clgeti(), impnld(), imaccess() +int fstati() +int tbpsta(), tbcigi() + +pointer imt, tnt # pointers for filename templates +int nin, nout # numbers of names in lists +pointer imtopen(), tbnopen() +int imtlen(), imtgetim(), tbnlen(), tbnget() + +begin + call smark (sp) + call salloc (inlist, SZ_LINE, TY_CHAR) + call salloc (outlist, SZ_LINE, TY_CHAR) + call salloc (temp, SZ_FNAME, TY_CHAR) + + # Get the names of the input tables. + if (fstati (STDIN, F_REDIR) == YES) + call strcpy ("STDIN", Memc[inlist], SZ_FNAME) + else + call clgstr ("intable", Memc[inlist], SZ_FNAME) + tnt = tbnopen (Memc[inlist]) + nin = tbnlen (tnt) + + # Get the names of the output images. + call clgstr ("output", Memc[outlist], SZ_FNAME) + imt = imtopen (Memc[outlist]) + nout = imtlen (imt) + + if (nin == 0) + call error (1, "no input table specified") + if (nout == 0) + call error (1, "no output image specified") + if (nin != nout) + call error (1, "input and output lists must have the same length") + + call clgstr ("colname", colname, SZ_COLNAME) + + # ndim is either zero or the dimension for new output images. + ndim = clgeti ("ndim") + if (ndim < 1) + ndim = 1 + do k = 1, IM_MAXDIM # initial values + axlen[k] = 1 + # Get the length of all but the last axis. + do k = 1, ndim-1 { + call sprintf (Memc[temp], SZ_FNAME, "n%d") + call pargi (k) + axlen[k] = clgeti (Memc[temp]) + } + + # Loop over the list of input tables. + while (tbnget (tnt, intable, SZ_FNAME) != EOF) { + + junk = imtgetim (imt, output, SZ_FNAME) + + tp = tbtopn (intable, READ_ONLY, NULL) + call tbcfnd (tp, colname, cp, 1) # only one column name + if (cp == NULL) { + call tbtclo (tp) + call error (1, "column not found") + } + nrows = tbpsta (tp, TBL_NROWS) + + # Open the output image. + if (imaccess (output, READ_WRITE) == YES) { + new_image = false + im = immap (output, READ_WRITE, NULL) + } else { + new_image = true + im = immap (output, NEW_IMAGE, NULL) + } + + if (new_image) { + # Set the size of the new image. + IM_NDIM(im) = ndim + npix = 1 # initial value + do k = 1, ndim-1 { + IM_LEN(im,k) = axlen[k] + npix = npix * axlen[k] + } + axlen[ndim] = nrows / npix + IM_LEN(im,ndim) = axlen[ndim] + + # The image data type is the same as that of the column. + dtype = tbcigi (cp, TBL_COL_DATATYPE) + if (dtype == TY_BOOL) + dtype = TY_SHORT + IM_PIXTYPE(im) = dtype + } + + nlines = 1 # initial value + do k = 2, IM_NDIM(im) + nlines = nlines * IM_LEN(im,k) + if (IM_LEN(im,1) * nlines != nrows) { + call tbtclo (tp) + call imunmap (im) + if (new_image) { + call imdelete (output) + call error (1, + "specified axis lengths are not consistent with size of table") + } else { + call error (1, + "size of existing image is not consistent with size of table") + } + } + + # Allocate space for the array of null flags (which we ignore), + # one element for each pixel in a line. + call salloc (nullflag, IM_LEN(im,1), TY_BOOL) + + # Copy the column into the image, one line at a time. + do k = 1, IM_MAXDIM + v[k] = 1 + frow = 1 + lrow = IM_LEN(im,1) + do k = 1, nlines { + junk = impnld (im, xp, v) + call tbcgtd (tp, cp, Memd[xp], Memb[nullflag], frow, lrow) + do i = 0, lrow-frow { + if (Memb[nullflag+i]) + Memd[xp+i] = -999.d0 + } + frow = frow + IM_LEN(im,1) + lrow = lrow + IM_LEN(im,1) + } + + call imunmap (im) + call tbtclo (tp) + } + + call sfree (sp) +end |