Initial commit

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