Initial commit

20 files changed, 4443 insertions, 0 deletions

diff --git a/pkg/obsolete/fits/README b/pkg/obsolete/fits/README
new file mode 100644
index 00000000..b9c10c7b
--- /dev/null
+++ b/pkg/obsolete/fits/README
@@ -0,0 +1 @@
+The IRAF  code for the FITS I/O package.
diff --git a/pkg/obsolete/fits/fits_cards.x b/pkg/obsolete/fits/fits_cards.x
new file mode 100644
index 00000000..fe4829c9
--- /dev/null
+++ b/pkg/obsolete/fits/fits_cards.x
@@ -0,0 +1,250 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include <imhdr.h>
+include "wfits.h"
+
+# WFT_STANDARD_CARD -- Procedure for fetching the minimum header
+# parameters required by fits. The end card is encoded separately.
+
+int procedure wft_standard_card (cardno, im, fits, axisno, card)
+
+int	cardno		# number of FITS standard card
+pointer	im		# pointer to the IRAF image
+pointer	fits		# pointer to the FITS structure
+int	axisno		# axis number
+char	card[ARB]	# FITS card image
+
+char	keyword[LEN_KEYWORD]
+errchk	wft_encodeb, wft_encodei, wft_encodel, wft_encode_axis
+include "wfits.com"
+
+begin
+	# Get mandatory keywords.
+	switch (cardno) {
+	case FIRST_CARD:
+	    call wft_encodeb ("SIMPLE", YES, card, "FITS STANDARD")
+	case SECOND_CARD:
+	    call wft_encodei ("BITPIX", FITS_BITPIX(fits), card,
+	    "FITS BITS/PIXEL")
+	case THIRD_CARD:
+	    call wft_encodei ("NAXIS", NAXIS(im), card, "NUMBER OF AXES")
+	default:
+	    call wft_encode_axis ("NAXIS", keyword, axisno)
+	    call wft_encodel (keyword, NAXISN(im, axisno), card, "")
+	    axisno = axisno + 1
+	}
+
+	return (YES)
+end
+
+
+# WFT_OPTION_CARD -- Procedure for fetching  optional FITS header parameters.
+# At present these are bscale, bzero, bunit, blank, object, origin, date,
+# irafmax, irafmin, iraf type and iraf bits per pixel. Blank is only encoded
+# if there are a nonzero number of blanks in the IRAF image. Bunit and object
+# are only encoded if the appropriate IRAF strings are defined. Bzero, bscale,
+# irafmax, irafmin, iraf type and iraf bits per pixel are only encoded if
+# there is a pixel file.
+
+int procedure wft_option_card (im, fits, optiono, card)
+
+pointer	im		# pointer to the IRAF image
+pointer	fits		# pointer to FITS structure
+int	optiono		# number of the option card
+char	card[ARB]	# FITS card image
+
+char	datestr[LEN_STRING]
+int	len_object, stat
+int	strlen()
+errchk	wft_encoded, wft_encodec, wft_encode_blank, wft_encoder, wft_encodei
+errchk	wft_encode_date 
+
+begin
+	 stat = YES
+
+	 # get optional keywords
+	 switch (optiono) {
+	 case KEY_BSCALE:
+	    if ((NAXIS(im) <= 0) || (FITS_BITPIX(fits) < 0))
+		stat = NO
+	    else {
+	        call wft_encoded ("BSCALE", BSCALE(fits), card,
+	        "REAL = TAPE*BSCALE + BZERO", NDEC_DOUBLE)
+	    }
+	 case KEY_BZERO:
+	    if ((NAXIS(im) <= 0) || (FITS_BITPIX(fits) < 0))
+		stat = NO
+	    else
+	        call wft_encoded ("BZERO", BZERO(fits), card, "", NDEC_DOUBLE)
+	 case KEY_BUNIT:
+	    stat = NO
+	 case KEY_BLANK:
+	    stat = NO
+	    #if (NBPIX(im) == 0)
+		#stat = NO
+	    #else
+		#call wft_encode_blank ("BLANK", BLANK_STRING(fits), card,
+		#"TAPE VALUE OF BLANK PIXEL")
+	case KEY_OBJECT:
+	    if (OBJECT(im) == EOS)
+		stat = NO
+	    else {
+		len_object = max (min (LEN_OBJECT, strlen (OBJECT(im))),
+		    LEN_STRING)
+		call wft_encodec ("OBJECT", OBJECT(im), len_object, card, "")
+	    }
+	case KEY_ORIGIN:
+	    call wft_encodec ("ORIGIN", "KPNO-IRAF", LEN_ORIGIN, card, "")
+	case KEY_DATE:
+	    call wft_encode_date (datestr, LEN_STRING)
+	    call wft_encodec ("DATE", datestr, LEN_STRING, card, "")
+	case KEY_IRAFNAME:
+	    len_object = max (min (LEN_OBJECT, strlen (IRAFNAME(fits))),
+		LEN_STRING)
+	    call wft_encodec ("IRAFNAME", IRAFNAME(fits), len_object, card, 
+		"NAME OF IRAF IMAGE FILE")
+	case KEY_IRAFMAX:
+	    if (NAXIS(im) <= 0)
+		stat = NO
+	    else
+	        call wft_encoder ("IRAF-MAX", IRAFMAX(fits), card, "DATA MAX",
+				NDEC_REAL)
+	case KEY_IRAFMIN:
+	    if (NAXIS(im) <= 0)
+		stat = NO
+	    else
+	        call wft_encoder ("IRAF-MIN", IRAFMIN(fits), card, "DATA MIN",
+		NDEC_REAL)
+	case KEY_IRAFBP:
+	    if (NAXIS(im) <= 0)
+		stat = NO
+	    else
+	        call wft_encodei ("IRAF-BPX", DATA_BITPIX(fits), card,
+	    	    "DATA BITS/PIXEL")
+	case KEY_IRAFTYPE:
+	    if (NAXIS(im) <= 0)
+		stat = NO
+	    else
+	        call wft_encodec ("IRAFTYPE", TYPE_STRING(fits), LEN_STRING,
+	    		    card, "PIXEL TYPE")
+	default:
+	    stat = NO
+	}
+
+	optiono = optiono + 1
+
+	return (stat)
+end
+
+
+# WFT_HISTORY_CARD -- Procedure to fetch a single history line, trim newlines
+# and pad with blanks to size LEN_CARD in order to create a FITS HISTORY card.
+
+int procedure wft_history_card (im, hp, card)
+
+pointer	im		# pointer to the IRAF image
+int	hp		# pointer to first character to extract from string
+char	card[ARB]	# FITS card image
+
+char	cval
+char	chfetch()
+
+begin
+	if (chfetch (HISTORY(im), hp, cval) == EOS)
+	    return (NO)
+	else {
+	    hp = hp - 1
+	    call strcpy ("HISTORY ", card, LEN_KEYWORD)
+	    call wft_fits_card (HISTORY(im), hp, card, COL_VALUE - 2, LEN_CARD,
+		'\n')
+	    return (YES)
+	}
+end
+
+
+# WFT_UNKNOWN_CARD  -- Procedure to fetch a single unknown
+# "line", trim newlines and pad blanks to size LEN_CARD in order to
+# create an unknown keyword card. At present user area information is
+# assumed to be in the form of FITS card images, less then or equal to
+# 80 characters and delimited by a newline.
+
+int procedure wft_unknown_card (im, up, card)
+
+pointer	im		# pointer to the IRAF image
+int	up		# pointer to next character in the unknown string
+char	card[ARB]	# FITS card image
+
+char	cval
+int	stat, axis, index
+char	chfetch()
+int	strmatch(), ctoi()
+
+begin
+	if (chfetch (UNKNOWN(im), up, cval) == EOS)
+	    return (NO)
+	else {
+	    up = up - 1
+	    stat = NO
+	    while (stat == NO)  {
+	        call wft_fits_card (UNKNOWN(im), up, card, 1, LEN_CARD, '\n')
+		if (card[1] == EOS)
+		    break
+	        if (strmatch (card, "^GROUPS  ") != 0) {
+		    stat = NO
+	        } else if (strmatch (card, "^SIMPLE  ") != 0) {
+		    stat = NO
+	        } else if (strmatch (card, "^BITPIX  ") != 0) {
+		    stat = NO
+	        } else if (strmatch (card, "^NAXIS   ") != 0) {
+		    stat = NO
+	        } else if (strmatch (card, "^NAXIS") != 0) {
+		    index = LEN_NAXIS_KYWRD + 1
+		    if (ctoi (card, index, axis) > 0)
+		        stat = NO
+		    else
+			stat = YES
+	        } else if (strmatch (card, "^GCOUNT  ") != 0) {
+		    stat = NO
+	        } else if (strmatch (card, "^PCOUNT  ") != 0) {
+		    stat = NO
+	        } else if (strmatch (card, "^PSIZE   ") != 0) {
+		    stat = NO
+	        } else if (strmatch (card, "^BSCALE  ") != 0) {
+		    stat = NO
+	        } else if (strmatch (card, "^BZERO   ") != 0) {
+		    stat = NO
+	        } else if (strmatch (card, "^BLANK   ") != 0) {
+		    stat = NO
+	        } else if (strmatch (card, "^IRAF-MAX") != 0) {
+		    stat = NO
+	        } else if (strmatch (card, "^IRAF-MIN") != 0) {
+		    stat = NO
+	        } else if (strmatch (card, "^IRAFTYPE") != 0) {
+		    stat = NO
+	        } else if (strmatch (card, "^IRAF-B/P") != 0) {
+		    stat = NO
+	        } else if (strmatch (card, "^IRAF-BPX") != 0) {
+		    stat = NO
+	        } else if (strmatch (card, "^END     ") != 0) {
+		    stat = NO
+	        } else
+                    stat = YES
+	    }
+
+	    return (stat)
+	}
+end
+
+
+# WFT_LAST_CARD -- Procedure to encode the FITS end card.
+
+int procedure wft_last_card (card)
+
+char	card[ARB]	# FITS card image
+
+begin
+	call sprintf (card, LEN_CARD, "%-8.8s  %70w")
+	    call pargstr ("END")
+
+	return (YES)
+end
diff --git a/pkg/obsolete/fits/fits_params.x b/pkg/obsolete/fits/fits_params.x
new file mode 100644
index 00000000..e3bf04ba
--- /dev/null
+++ b/pkg/obsolete/fits/fits_params.x
@@ -0,0 +1,234 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include <time.h>
+include "wfits.h"
+
+# WFT_ENCODEB -- Procedure to encode a boolean parameter into a FITS card.
+
+procedure wft_encodeb (keyword, param, card, comment)
+
+char	keyword[ARB]	# FITS keyword
+int	param		# integer parameter equal to YES/NO
+char	card[ARB]	# FITS card image
+char	comment[ARB]	# FITS comment string
+
+char	truth
+
+begin
+	if (param == YES)
+	    truth = 'T'
+	else
+	    truth = 'F'
+
+	call sprintf (card, LEN_CARD, "%-8.8s= %20c  /  %-45.45s")
+	    call pargstr (keyword)
+	    call pargc (truth)
+	    call pargstr (comment)
+end
+
+
+# WFT_ENCODEI -- Procedure to encode an integer parameter into a FITS card.
+
+procedure wft_encodei (keyword, param, card, comment)
+
+char	keyword[ARB]	# FITS keyword
+int	param		# integer parameter
+char	card[ARB]	# FITS card image
+char	comment[ARB]	# FITS comment string
+
+begin
+	call sprintf (card, LEN_CARD, "%-8.8s= %20d  /  %-45.45s")
+	    call pargstr (keyword)
+	    call pargi (param)
+	    call pargstr (comment)
+end
+
+
+# WFT_ENCODEL -- Procedure to encode a long parameter into a FITS card.
+
+procedure wft_encodel (keyword, param, card, comment)
+
+char	keyword[ARB]		# FITS keyword
+long	param			# long integer parameter
+char	card[ARB]		# FITS card image
+char	comment[ARB]		# FITS comment string
+
+begin
+	call sprintf (card, LEN_CARD, "%-8.8s= %20d  /  %-45.45s")
+	    call pargstr (keyword)
+	    call pargl (param)
+	    call pargstr (comment)
+end
+
+
+# WFT_ENCODER -- Procedure to encode a real parameter into a FITS card.
+
+procedure wft_encoder (keyword, param, card, comment, precision)
+
+char	keyword[ARB]		# FITS keyword
+real	param			# real parameter
+char	card[ARB]		# FITS card image
+char	comment[ARB]		# FITS comment card
+int	precision		# precision of real
+
+begin
+	call sprintf (card, LEN_CARD, "%-8.8s= %20.*e  /  %-45.45s")
+	    call pargstr (keyword)
+	    call pargi (precision)
+	    call pargr (param)
+	    call pargstr (comment)
+end
+
+
+# WFT_ENCODED -- Procedure to encode a double parameter into a FITS card.
+
+procedure wft_encoded (keyword, param, card, comment, precision)
+
+char	keyword[ARB]		# FITS keyword
+double	param			# double parameter
+char	card[ARB]		# FITS card image
+char	comment[ARB]		# FITS comment string
+int	precision		# FITS precision
+
+begin
+	call sprintf (card, LEN_CARD, "%-8.8s= %20.*e  /  %-45.45s")
+	    call pargstr (keyword)
+	    call pargi (precision)
+	    call pargd (param)
+	    call pargstr (comment)
+end
+
+
+# WFT_ENCODE_AXIS -- Procedure to add the axis number to axis dependent
+# keywords.
+
+procedure wft_encode_axis (root, keyword, axisno)
+
+char	root[ARB]		# FITS root keyword
+char	keyword[ARB]		# FITS keyword
+int	axisno			# FITS axis number
+
+begin
+	call strcpy (root, keyword, LEN_KEYWORD)
+	call sprintf (keyword, LEN_KEYWORD, "%-5.5s%-3.3s")
+	    call pargstr (root)
+	    call pargi (axisno)
+end
+
+
+# WFT_ENCODEC -- Procedure to encode an IRAF string parameter into a FITS card.
+
+procedure wft_encodec (keyword, param, maxch, card, comment)
+
+char	keyword[ARB]	# FITS keyword
+char	param[ARB]	# FITS string parameter
+int	maxch		# maximum number of characters in string parameter
+char	card[ARB]	# FITS card image
+char	comment[ARB]	# comment string
+
+char	strparam[LEN_ALIGN+2]
+int	maxchar, nblanks
+
+begin
+	maxchar = min (maxch, LEN_OBJECT)
+	if (maxchar <= LEN_ALIGN - 1) {
+	    strparam[1] = '\''
+	    call sprintf (strparam[2], maxchar, "%*.*s")
+		call pargi (-maxchar)
+		call pargi (maxchar)
+		call pargstr (param)
+	    strparam[maxchar+2] = '\''
+	    strparam[maxchar+3] = EOS
+            call sprintf (card, LEN_CARD, "%-8.8s= %-20.20s  /  %-45.45s")
+	        call pargstr (keyword)
+	        call pargstr (strparam)
+	        call pargstr (comment)
+	} else {
+	    nblanks = LEN_OBJECT - maxchar
+            call sprintf (card, LEN_CARD, "%-8.8s= '%*.*s'  /  %*.*s")
+	        call pargstr (keyword)
+	        call pargi (-maxchar)
+	        call pargi (maxchar)
+	        call pargstr (param)
+	        call pargi (-nblanks)
+	        call pargi (nblanks)
+	        call pargstr (comment)
+	}
+end
+
+
+# WFT_ENCODE_BLANK -- Procedure to encode the FITS blank parameter. Necessary
+# because the 32 bit blank value equals INDEFL.
+
+procedure wft_encode_blank (keyword, blank_str, card, comment)
+
+char	keyword[ARB]		# FITS keyword
+char	blank_str[ARB]		# string containing values of FITS blank integer
+char	card[ARB]		# FITS card image
+char	comment[ARB]		# FITS comment string
+
+begin
+    call sprintf (card, LEN_CARD, "%-8.8s= %20.20s  /  %-45.45s")
+	call pargstr (keyword)
+	call pargstr (blank_str)
+	call pargstr (comment)
+end
+
+
+# WFT_ENCODE_DATE -- Procedure to encode the date in the form dd-mm-yy.
+
+procedure wft_encode_date (datestr, szdate)
+
+char	datestr[ARB]	# string containing the date
+int	szdate		# number of chars in the date string
+
+long	ctime
+int	time[LEN_TMSTRUCT]
+long	clktime()
+
+begin
+	ctime = clktime (long (0))
+	call brktime (ctime, time)
+
+	call sprintf (datestr, szdate, "%02s-%02s-%02s")
+	    call pargi (TM_MDAY(time))
+	    call pargi (TM_MONTH(time))
+	    call pargi (mod (TM_YEAR(time), CENTURY))
+end
+
+
+# WFT_FITS_CARD --  Procedure to fetch a single line from a string parameter
+# padding it to a maximum of maxcols characters and trimmimg the delim
+# character.
+
+procedure wft_fits_card (instr, ip, card, col_out, maxcols, delim)
+
+char	instr[ARB]	# input string
+int	ip		# input string pointer, updated at each call
+char	card[ARB]	# FITS card image
+int	col_out		# pointer to column in card
+int	maxcols		# maximum columns in card
+int	delim		# 1 character string delimiter
+
+int	op
+
+begin
+	op = col_out
+
+	# Copy string
+	while (op <= maxcols && instr[ip] != EOS && instr[ip] != delim) {
+	    card[op] = instr[ip]
+	    ip = ip + 1
+	    op = op + 1
+	}
+
+	# Fill remainder of card with blanks
+	while (op <= maxcols ) {
+	    card[op] = ' '
+	    op = op + 1
+	}
+
+	if (instr[ip] == delim)
+	    ip = ip + 1
+
+end
diff --git a/pkg/obsolete/fits/fits_read.x b/pkg/obsolete/fits/fits_read.x
new file mode 100644
index 00000000..bcbdb745
--- /dev/null
+++ b/pkg/obsolete/fits/fits_read.x
@@ -0,0 +1,173 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include <error.h>
+include <imhdr.h>
+include <fset.h>
+include	"rfits.h"
+
+# RFT_READ_FITZ -- Convert a FITS file. An EOT is signalled by returning EOF.
+
+int procedure rft_read_fitz (fitsfile, iraffile)
+
+char	fitsfile[ARB]		# FITS file name
+char	iraffile[ARB]		# IRAF file name
+
+int	fits_fd, stat, min_lenuserarea, ip
+pointer	im, sp, fits, envstr
+int	rft_read_header(), mtopen(), immap(), strlen(), envfind(), ctoi()
+errchk	smark, sfree, salloc, rft_read_header, rft_read_image, rft_find_eof()
+errchk	rft_scan_file, mtopen, immap, imdelete, close, imunmap
+
+include	"rfits.com"
+
+begin
+	# Open input FITS data.
+	fits_fd = mtopen (fitsfile, READ_ONLY, 0)
+
+	# Allocate memory for program data structure.
+	call smark (sp)
+	call salloc (fits, LEN_FITS, TY_STRUCT)
+	call salloc (envstr, SZ_FNAME, TY_CHAR)
+
+	# Set up for printing a long or a short header.
+	if (long_header == YES || short_header == YES) {
+	    if (make_image == YES) {
+	        call printf ("File: %s  ")
+		    call pargstr (iraffile)
+	    } else {
+		call printf ("File: %s  ")
+		    call pargstr (fitsfile)
+	    }
+	    if (long_header == YES)
+		call printf ("\n")
+	}
+	call flush (STDOUT)
+
+	# Create the IRAF image header. If only a header listing is desired
+	# then map the scratch image onto DEV$NULL (faster than a real file).
+
+	if (make_image == NO)
+	    call strcpy ("dev$null", iraffile, SZ_FNAME)
+	if (envfind ("min_lenuserarea", Memc[envstr], SZ_FNAME) > 0) {
+	    ip = 1
+	    if (ctoi (Memc[envstr], ip, min_lenuserarea) <= 0)
+		min_lenuserarea = LEN_USERAREA
+	    else
+		min_lenuserarea = max (LEN_USERAREA, min_lenuserarea)
+	} else
+	    min_lenuserarea = LEN_USERAREA
+	im = immap (iraffile, NEW_IMAGE, min_lenuserarea)
+
+	# Read header.  EOT is signalled by an EOF status from fits_read_header.
+	# Create an IRAF image if desired.
+
+	iferr {
+	    IRAFNAME(fits) = EOS
+	    stat = rft_read_header (fits_fd, fits, im)
+	    if (stat == EOF)
+	        call printf ("End of data\n")
+	    else {
+	        if (make_image == YES) {
+	            call rft_read_image (fits_fd, fits, im)
+		    if (fe > 0.0)
+			call rft_find_eof (fits_fd)
+		} else if (fe > 0.0)
+		    call rft_scan_file (fits_fd, fits, im, fe)
+	    }
+	} then {
+	    call flush (STDOUT)
+	    call erract (EA_WARN)
+	}
+
+	# Close files and clean up.
+	call imunmap (im)
+
+	# Optionally restore the old IRAF name.
+	if (stat == EOF || make_image == NO) {
+	    call imdelete (iraffile)
+	} else if (old_name == YES && strlen (IRAFNAME(fits)) != 0) {
+	    iferr {
+		call imgimage (IRAFNAME(fits), IRAFNAME(fits), SZ_FNAME)
+	        call imrename (iraffile, IRAFNAME(fits))
+	    } then {
+		call printf ("   Cannot rename image %s to %s\n")
+		    call pargstr (iraffile)
+		    call pargstr (IRAFNAME(fits))
+	        call flush (STDOUT)
+	        call erract (EA_WARN)
+	    } else {
+	        call printf ("    File: %s restored to IRAF File: %s\n")
+		    call pargstr (iraffile)
+		    call pargstr (IRAFNAME(fits))
+	    }
+	}
+
+	if (long_header == YES)
+	    call printf ("\n")
+
+	call close (fits_fd)
+	call sfree (sp)
+
+	return (stat)
+end
+
+
+# RFT_FIND_EOF -- Read the FITS data file until EOF is reached.
+
+procedure rft_find_eof (fd)
+
+int	fd			# the FITS file descriptor
+
+int	szbuf
+pointer	sp, buf
+int	fstati(), read()
+errchk	read
+
+begin
+	# Scan through the file.
+	szbuf = fstati (fd, F_BUFSIZE)
+	call smark (sp)
+	call salloc (buf, szbuf, TY_CHAR)
+	while (read (fd, Memc[buf], szbuf) != EOF)
+	    ;
+	call sfree (sp)
+end
+
+
+# RFT_SCAN_FILE -- Determine whether it is more efficient to read the
+# entire file or to skip forward to the next file if the parameter
+# make_image was set to no.
+
+procedure  rft_scan_file (fd, fits, im, fe)
+
+int	fd			# the FITS file descriptor
+pointer	fits			# pointer to the FITS descriptor
+pointer	im			# pointer to the output image
+real	fe			# maximum file size in Kb for scan mode
+
+int	i, szbuf
+pointer	sp, buf
+real	file_size
+int	fstati(), read()
+errchk	read
+
+begin
+	# Compute the file size in Kb and return if it is bigger than fe.
+	file_size = 1.0
+	do i = 1, IM_NDIM(im)
+	    file_size = file_size * IM_LEN(im,i)
+	if (IM_NDIM(im) <= 0)
+	    file_size = 0.0
+	else
+	    file_size = file_size * abs (BITPIX(fits)) / FITS_BYTE / 1.0e3
+	if (file_size >= fe)
+	    return
+
+	# Scan through the file.
+	szbuf = fstati (fd, F_BUFSIZE)
+	call smark (sp)
+	call salloc (buf, szbuf, TY_CHAR)
+	while (read (fd, Memc[buf], szbuf) != EOF)
+	    ;
+	call sfree (sp)
+end
diff --git a/pkg/obsolete/fits/fits_rheader.x b/pkg/obsolete/fits/fits_rheader.x
new file mode 100644
index 00000000..8efb5562
--- /dev/null
+++ b/pkg/obsolete/fits/fits_rheader.x
@@ -0,0 +1,575 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include <ctype.h>
+include <imhdr.h>
+include <imio.h>
+include <mach.h>
+include	"rfits.h"
+
+define	NEPSILON 10.0d0    # number of machine epsilon
+
+# RFT_READ_HEADER -- Read a FITS header.
+# If BSCALE and BZERO are different from 1.0 and 0.0  scale is set to true
+# otherwise scale is false.
+# EOT is detected by an EOF on the first read and EOF is returned to the calling
+# routine.  Any error is passed to the calling routine.
+
+int procedure rft_read_header (fits_fd, fits, im)
+
+int	fits_fd			# FITS file descriptor
+pointer	fits			# FITS data structure
+pointer	im			# IRAF image descriptor
+
+int	i, stat, nread, max_lenuser, fd_usr, ndiscard
+char	card[LEN_CARD+1]
+int	rft_decode_card(), rft_init_read_pixels(), rft_read_pixels(), strmatch()
+int	stropen()
+errchk	rft_decode_card, rft_init_read_pixels, rft_read_pixels
+errchk	stropen, close
+
+include "rfits.com"
+
+begin
+	# Initialization.
+	SIMPLE(fits) = NO
+	BITPIX(fits) = INDEFI
+	NAXIS(im) = INDEFI
+	do i = 1, IM_MAXDIM
+	    IM_LEN(im,i) = INDEFL
+	SCALE(fits) = NO
+	FITS_BSCALE(fits) = 1.0d0
+	FITS_BZERO(fits) = 0.0d0
+	BLANKS(fits) = NO
+	BLANK_VALUE(fits) = INDEFL
+	NRECORDS(fits) = 0
+	ndiscard = 0
+	max_lenuser = (LEN_IMDES + IM_LENHDRMEM(im) - IMU) * SZ_STRUCT - 1
+
+	# The FITS header is character data in FITS_BYTE form. Open the
+	# header for reading. Open the user area which is a character
+	# string as a file.
+
+	i = rft_init_read_pixels (len_record, FITS_BYTE, LSBF, TY_CHAR)
+	fd_usr = stropen (UNKNOWN(im), max_lenuser, NEW_FILE)
+
+	# Loop until the END card is encountered.
+	nread = 0
+	repeat {
+
+	    # Read the card.
+	    i = rft_read_pixels (fits_fd, card, LEN_CARD, NRECORDS(fits), 1)
+	    card[LEN_CARD + 1] = '\n'
+	    card[LEN_CARD + 2] = EOS
+
+	    # Decode the card images.
+	    if ((i == EOF) && (nread == 0)) {
+		return (EOF)
+	    } else if ((nread == 0) && strmatch (card, "^SIMPLE  ") == 0) {
+		call flush (STDOUT)
+		call error (30,
+		    "RFT_READ_HEADER: Not a FITS file (no SIMPLE keyword)")
+	    } else if (i != LEN_CARD) {
+	        call error (2, "RFT_READ_HEADER: Error reading FITS header")
+	    } else
+	        nread = nread + 1
+
+	    # Remove contaminating control characters and replace with blanks.
+	    call rft_control_to_blank (card, card, LEN_CARD)
+
+	    # Print FITS card images if long_header option specified.
+	    if (long_header == YES) {
+		call printf ("%-80.80s\n")
+		    call pargstr (card)
+	    }
+
+	    # Stat = YES if FITS END card is encountered.
+	    stat = rft_decode_card (fits, im, fd_usr, card, ndiscard)
+
+	} until (stat == YES)
+
+	# Print optional short header.
+	if (short_header == YES && long_header == NO) {
+	    if (make_image == NO && old_name == YES) {
+	        call printf ("(%s)  %-20.20s  ")
+		    call pargstr (IRAFNAME(fits))
+		    call pargstr (OBJECT(im))
+	    } else {
+	        call printf ("%-20.20s  ")
+		    call pargstr (OBJECT(im))
+	    }
+	    do i = 1, NAXIS(im) {
+		if (i == 1) {
+		    call printf ("Size = %d")
+		        call pargl (NAXISN(im,i))
+		} else {
+		    call printf (" x %d")
+		        call pargl (NAXISN(im,i))
+		}
+	    }
+	    call printf ("\n")
+	}
+
+	# Let the user know if there is not enough space in the user area.
+	if (ndiscard > 0) {
+	    call printf (
+	        "Warning: User area too small %d card images discarded\n")
+	    call pargi (ndiscard)
+	    call rft_last_user (UNKNOWN(im), max_lenuser)
+	}
+
+	call close (fd_usr)
+	return (OK)
+end
+
+
+# RFT_CONTROL_TO_BLANK -- Replace an ACSII control characters in the
+# FITS card image with blanks.
+
+procedure rft_control_to_blank (incard, outcard, len_card)
+
+char	incard[ARB]		# the input FITS card image
+char	outcard[ARB]		# the output FITS card image
+int	len_card		# the length of the FITS card image
+
+int	i
+
+begin
+	for (i = 1; i <= len_card; i = i + 1) {
+	    if (IS_PRINT(incard[i]))
+	        outcard[i] = incard[i]
+	    else
+		outcard[i] = ' '
+	}
+end
+
+
+# RFT_DECODE_CARD -- Decode a FITS card and return YES when the END
+# card is encountered.  The keywords understood are given in rfits.h.
+
+int procedure rft_decode_card (fits, im, fd_usr, card, ndiscard)
+
+pointer	fits		# FITS data structure
+pointer	im		# IRAF image descriptor
+int	fd_usr		# file descriptor of user area
+char	card[ARB]	# FITS card
+int	ndiscard	# Number of cards for which no space available
+
+char	cval, str[LEN_CARD]
+double	dval
+int	nchar, ival, i, j, k, len
+pointer	sp, comment
+
+bool	rft_equald()
+int	strmatch(), ctoi(), ctol(), ctod(), cctoc(), rft_hms()
+errchk	putline
+
+include	"rfits.com"
+
+begin
+	call smark (sp)
+	call salloc (comment, SZ_LINE, TY_CHAR)
+
+	i = COL_VALUE
+	if (strmatch (card, "^END     ") != 0) {
+	    call sfree (sp)
+	    return(YES)
+	} else if (strmatch (card, "^SIMPLE  ") != 0) {
+	    if (SIMPLE(fits) == YES)
+		call printf ("Warning: Duplicate SIMPLE keyword ignored\n")
+	    else {
+	        nchar = cctoc (card, i, cval)
+	        if (cval != 'T')
+		    call error (13, "RFT_DECODE_CARD: Non-standard FITS format")
+		SIMPLE(fits) = YES
+	    }
+	} else if (strmatch (card, "^BITPIX  ") != 0) {
+	    if (! IS_INDEFI(BITPIX(fits)))
+		call printf ("Warning: Duplicate BITPIX keyword ignored\n")
+	    else
+	        nchar = ctoi (card, i, BITPIX(fits))
+	} else if (strmatch (card, "^BLANK   ") != 0) {
+	    BLANKS(fits) = YES
+	    nchar = ctol (card, i, BLANK_VALUE(fits))
+	} else if (strmatch (card, "^NAXIS   ") != 0) {
+	    if (! IS_INDEFI(NAXIS(im)))
+		call printf ("Warning: Duplicate NAXIS keyword ignored\n")
+	    else
+	        nchar = ctoi (card, i, NAXIS(im))
+	    if (NAXIS(im) > IM_MAXDIM)
+		call error (5, "RFT_DECODE_CARD: FITS NAXIS too large")
+	} else if (strmatch (card, "^NAXIS") != 0) {
+	    k = strmatch (card, "^NAXIS")
+	    nchar = ctoi (card, k, j)
+	    if (! IS_INDEFL(NAXISN(im,j))) {
+		call printf ("Warning: Duplicate NAXIS%d keyword ignored\n")
+		    call pargi (j)
+	    } else
+	        nchar = ctol (card, i, NAXISN(im, j))
+	} else if (strmatch (card, "^GROUPS  ") != 0) {
+	    nchar = cctoc (card, i, cval)
+	    if (cval == 'T') {
+		NAXIS(im) = 0
+		call error (6, "RFT_DECODE_CARD: Group data not implemented")
+	    }
+	} else if (strmatch (card, "^TABLES  ") != 0) {
+	    nchar = ctoi (card, i, ival)
+	    if (ival > 0)
+		call printf ("Warning: FITS special records not decoded\n")
+	} else if (strmatch (card, "^BSCALE  ") != 0) {
+	    nchar = ctod (card, i, dval)
+	    if (nchar > 0)
+	        FITS_BSCALE(fits) = dval
+	    else
+		call printf ("Warning: Error decoding BSCALE, BSCALE=1.0\n")
+	    if (! rft_equald (dval, 1.0d0) && (scale == YES))
+		SCALE(fits) = YES
+	} else if (strmatch (card, "^BZERO   ") != 0) {
+	    nchar = ctod (card, i, dval)
+	    if (nchar > 0)
+	        FITS_BZERO(fits) = dval
+	    else
+		call printf ("Warning: Error decoding BZERO, BZERO=0.0\n")
+	    if (! rft_equald (dval, 0.0d0) && (scale == YES))
+		SCALE(fits) = YES
+	} else if (strmatch (card, "^OBJECT  ") != 0) {
+	    call rft_get_fits_string (card, OBJECT(im), SZ_OBJECT)
+	} else if (strmatch (card, "^IRAFNAME") != 0) {
+	    call rft_get_fits_string (card, IRAFNAME(fits), SZ_FNAME)
+	} else if (strmatch (card, "^ORIGIN  ") != 0) {
+	    call rft_trim_card (card, card, LEN_CARD)
+	    call strcat (card[i], HISTORY(im), SZ_HISTORY)
+	} else if (strmatch (card, "^DATE    ") !=  0) {
+	    call rft_trim_card (card, card, LEN_CARD)
+	    call strcat (card[i], HISTORY(im), SZ_HISTORY)
+	#} else if (strmatch (card, "^HISTORY ") != 0) {
+	    #call rft_trim_card (card, card, LEN_CARD)
+	    #call strcat (card[i - 2], HISTORY(im), SZ_HISTORY)
+	} else if (strmatch (card, "^UT      ") != 0) {
+	    len = rft_hms (card, str, Memc[comment], LEN_CARD)
+	    if (len > 0) {
+		call wft_encodec ("UT", str, len, card, Memc[comment])
+		card[LEN_CARD+1] = '\n'
+		card[LEN_CARD+2] = EOS
+	    }
+	    if (ndiscard > 1)
+		ndiscard = ndiscard + 1
+	    else {
+	        iferr (call putline (fd_usr, card))
+		    ndiscard = ndiscard + 1
+	    }
+	} else if (strmatch (card, "^ZD      ") != 0) {
+	    len = rft_hms (card, str, Memc[comment], LEN_CARD)
+	    if (len > 0) {
+		call wft_encodec ("ZD", str, len, card, Memc[comment])
+		card[LEN_CARD+1] = '\n'
+		card[LEN_CARD+2] = EOS
+	    }
+	    if (ndiscard > 1)
+		ndiscard = ndiscard + 1
+	    else {
+	        iferr (call putline (fd_usr, card))
+		    ndiscard = ndiscard + 1
+	    }
+	} else if (strmatch (card, "^ST      ") != 0) {
+	    len = rft_hms (card, str, Memc[comment], LEN_CARD)
+	    if (len > 0) {
+		call wft_encodec ("ST", str, len, card, Memc[comment])
+		card[LEN_CARD+1] = '\n'
+		card[LEN_CARD+2] = EOS
+	    }
+	    if (ndiscard > 1)
+		ndiscard = ndiscard + 1
+	    else {
+	        iferr (call putline (fd_usr, card))
+		    ndiscard = ndiscard + 1
+	    }
+	} else if (strmatch (card, "^RA      ") != 0) {
+	    len = rft_hms (card, str, Memc[comment], LEN_CARD)
+	    if (len > 0) {
+		call wft_encodec ("RA", str, len, card, Memc[comment])
+		card[LEN_CARD+1] = '\n'
+		card[LEN_CARD+2] = EOS
+	    }
+	    if (ndiscard > 1)
+		ndiscard = ndiscard + 1
+	    else {
+	        iferr (call putline (fd_usr, card))
+		    ndiscard = ndiscard + 1
+	    }
+	} else if (strmatch (card, "^DEC     ") != 0) {
+	    len = rft_hms (card, str, Memc[comment], LEN_CARD)
+	    if (len > 0) {
+		call wft_encodec ("DEC", str, len, card, Memc[comment])
+		card[LEN_CARD+1] = '\n'
+		card[LEN_CARD+2] = EOS
+	    }
+	    if (ndiscard > 1)
+		ndiscard = ndiscard + 1
+	    else {
+	        iferr (call putline (fd_usr, card))
+		    ndiscard = ndiscard + 1
+	    }
+	} else {
+	    if (ndiscard > 1)
+		ndiscard = ndiscard + 1
+	    else {
+	        iferr (call putline (fd_usr, card))
+		    ndiscard = ndiscard + 1
+	    }
+	}
+
+	call sfree (sp)
+
+	return (NO)
+
+end
+
+
+# RFT_HMS -- Procedure to decode a FITS HMS card from the mountain.
+
+int procedure rft_hms (card, str, comment, maxch)
+
+char	card[ARB]		# FITS card
+char	str[ARB]		# string
+char	comment[ARB]		# comment string
+int	maxch			# maximum number of characters
+
+char	colon, minus
+int	ip, nchar, fst, lst, deg, min
+real	sec
+int	stridx(), strldx(), strlen(), ctoi(), ctor()
+
+begin
+	# Return if not a FITS string parameter.
+	if (card[COL_VALUE] != '\'')
+	    return (0)
+
+	# Set up key characters.
+	colon = ':'
+	minus = '-'
+
+	# Get the FITS string.
+	call rft_get_fits_string (card, str, maxch)
+
+	# Get the comment string.
+	call rft_get_comment (card, comment, maxch)
+
+	# Test for blank string and for 2 colon delimiters.
+	if (str[1] == EOS)
+	    return (0)
+	fst = stridx (colon, str)
+	if (fst == 0)
+	    return (0)
+	lst = strldx (colon, str)
+	if (lst == 0)
+	    return (0)
+	if (fst == lst)
+	    return (0)
+
+	# Decode the degrees field.
+	ip = 1
+	while (IS_WHITE(str[ip]))
+	    ip = ip + 1
+	if (str[ip] == '+' || str[ip] == '-')
+	    ip = ip + 1
+	nchar = ctoi (str, ip, deg)
+	if (nchar == 0)
+	    deg = 0
+
+	# Decode the minutes field.
+	ip = fst + 1
+	while (IS_WHITE(str[ip]))
+	    ip = ip + 1
+	if (str[ip] == '+' || str[ip] == '-')
+	    ip = ip + 1
+	nchar = ctoi (str, ip, min)
+	if (nchar == 0)
+	    min = 0
+
+	# Decode the seconds field.
+	ip = lst + 1
+	while (IS_WHITE(str[ip]))
+	    ip = ip + 1
+	if (str[ip] == '+' || str[ip] == '-')
+	    ip = ip + 1
+	nchar = ctor (str, ip, sec)
+	if (nchar == 0)
+	    sec = 0.0
+
+	# Reformat the HMS card.
+	if (stridx (minus, str) > 0 || deg < 0 || min < 0 || sec < 0.0) {
+	    call sprintf (str, maxch, "%c%d:%02d:%05.2f")
+		call pargc (minus)
+	        call pargi (abs (deg))
+	        call pargi (abs (min))
+	        call pargr (abs (sec))
+	} else {
+	    call sprintf (str, maxch, "%2d:%02d:%05.2f")
+	        call pargi (deg)
+	        call pargi (abs (min))
+	        call pargr (abs (sec))
+	}
+
+	return (strlen (str))
+end
+
+
+# RFT_GET_COMMENT -- Extract the comment field from a FITS card.
+
+procedure rft_get_comment (card, comment, maxch)
+
+char	card[ARB]		# FITS card
+char	comment[ARB]		# comment string
+int	maxch			# maximum number of characters
+
+int	istart, j
+
+begin
+	istart = 0
+	for (j = LEN_CARD; (j >= 1) && (card[j] != '\''); j = j - 1) {
+	    if (card[j] == '/') {
+		for (istart = j + 1; IS_WHITE(card[istart]) && istart <=
+		    LEN_CARD; istart = istart + 1)
+		    ;
+		break
+	    }
+	}
+
+	if (istart == 0)
+	    comment[1] = EOS
+	else
+	    call strcpy (card[istart], comment, LEN_CARD - istart + 1 )
+end
+
+
+# RFT_GET_FITS_STRING -- Extract a string from a FITS card and trim trailing
+# blanks. The EOS is marked by either ', /, or the end of the card.
+# There may be an optional opening ' (FITS standard).
+
+procedure rft_get_fits_string (card, str, maxchar)
+
+char	card[ARB]		# FITS card
+char	str[ARB]		# FITS string
+int	maxchar			# maximum number of characters
+
+int	j, istart, nchar
+
+begin
+	# Check for opening quote
+	for (istart = COL_VALUE; istart <= LEN_CARD && card[istart] != '\'';
+	    istart = istart + 1)
+	    ;
+	istart = istart + 1
+
+	# Check for closing quote.
+	for (j = istart; (j<LEN_CARD)&&(card[j]!='\''); j = j + 1)
+	    ;
+	for (j = j - 1; (j >= istart) && (card[j] == ' '); j = j - 1)
+	    ;
+	nchar = min (maxchar, j - istart + 1)
+
+	# Copy the string.
+	if (nchar <= 0)
+	    str[1] = EOS
+	else
+	    call strcpy (card[istart], str, nchar)
+end
+
+
+# RFT_EQUALD -- Procedure to compare two double precision numbers for equality
+# to within the machine precision for doubles.
+
+bool procedure rft_equald (x, y)
+
+double	x, y		# the two numbers to be compared for equality
+
+int	ex, ey
+double	x1, x2, normed_x, normed_y
+
+begin
+	if (x == y)
+	    return (true)
+
+	call rft_normd (x, normed_x, ex)
+	call rft_normd (y, normed_y, ey)
+
+	if (ex != ey)
+	    return (false)
+	else {
+	    x1 = 1.0d0 + abs (normed_x - normed_y)
+	    x2 = 1.0d0 + NEPSILON * EPSILOND
+	    return (x1 <= x2)
+	}
+end
+
+
+# RFT_NORMED -- Normalize a double precision number x to the value normed_x,
+# in the range [1-10]. Expon is returned such that x = normed_x *
+# (10.0d0 ** expon).
+
+procedure rft_normd (x, normed_x, expon)
+
+double	x			# number to be normailized
+double	normed_x		# normalized number
+int	expon			# exponent
+
+double	ax
+
+begin
+	ax = abs (x)
+	expon = 0
+
+	if (ax > 0) {
+	    while (ax < (1.0d0 - NEPSILON * EPSILOND)) {
+		ax = ax * 10.0d0
+		expon = expon - 1
+	    }
+
+	    while (ax >= (10.0d0 - NEPSILON * EPSILOND)) {
+		ax = ax / 10.0d0
+		expon = expon + 1
+	    }
+	}
+
+	if (x < 0)
+	    normed_x = -ax
+	else
+	    normed_x = ax
+end
+
+
+# RFT_TRIM_CARD -- Procedure to trim trailing whitespace from the card
+
+procedure rft_trim_card (incard, outcard, maxch)
+
+char	incard[ARB]		# input FITS card image 
+char	outcard[ARB]		# output FITS card
+int	maxch			# maximum size of card
+
+int	ip
+
+begin
+	ip = maxch
+	while (incard[ip] == ' ' || incard[ip] == '\t' || incard[ip] == '\0')
+	    ip = ip - 1
+	call amovc (incard, outcard, ip)
+	outcard[ip+1] = '\n'
+	outcard[ip+2] = EOS
+end
+
+
+# RFT_LAST_CARD -- Remove a partially written card from the data base
+
+procedure rft_last_user (user, maxch) 
+
+char	user[ARB]	# user area
+int	maxch		# maximum number of characters
+
+int	ip
+
+begin
+	ip = maxch
+	while (user[ip] != '\n')
+	    ip = ip - 1
+	user[ip+1] = EOS
+end
diff --git a/pkg/obsolete/fits/fits_rimage.x b/pkg/obsolete/fits/fits_rimage.x
new file mode 100644
index 00000000..e74e2dbb
--- /dev/null
+++ b/pkg/obsolete/fits/fits_rimage.x
@@ -0,0 +1,605 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<imhdr.h>
+include <mach.h>
+include <fset.h>
+include	"rfits.h"
+
+# RFT_READ_IMAGE --  Convert FITS image pixels to IRAF image pixels.
+
+procedure rft_read_image (fits_fd, fits, im)
+
+int	fits_fd		# FITS file descriptor
+pointer	fits		# FITS data structure
+pointer	im		# IRAF image descriptor
+
+int	i, npix, npix_record, blksize, ndummy
+long	v[IM_MAXDIM], nlines, il
+pointer	tempbuf, buf
+real	linemax, linemin, lirafmin, lirafmax
+double	dblank
+
+long	clktime()
+int	fstati(), rft_init_read_pixels(), rft_read_pixels()
+
+errchk	malloc, mfree, rft_init_read_pixels, rft_read_pixels, rft_lscale_pix
+errchk	rft_lchange_pix, rft_rchange_pix, rfit_dchange_pix, rft_put_image_line
+errchk	rft_pix_limits, rft_rscale_pix, rft_dscale_pix
+
+include	"rfits.com"
+
+begin
+	# No pixel file was created.
+	if (NAXIS(im) == 0) {
+	    call printf ("Warning: No pixel file created\n")
+	    return
+	}
+
+	# Initialize the header.
+	call rft_set_image_header (fits, im)
+
+	# Compute the number of columns and lines in the image.
+	npix = NAXISN(im, 1)
+	nlines = 1
+	do i = 2, NAXIS(im)
+	    nlines = nlines * NAXISN(im, i)
+	lirafmax = -MAX_REAL
+	lirafmin = MAX_REAL
+
+	# Compute the number of pixels per record and the number of records
+	# per output block.
+
+	npix_record = len_record * FITS_BYTE / abs (BITPIX(fits))
+	blksize = fstati (fits_fd, F_SZBBLK)
+	if (mod (blksize, FITS_RECORD) == 0)
+	    blksize = blksize / FITS_RECORD
+	else
+	    blksize = 1
+
+	# FITS data is converted to type  LONG, REAL or DOUBLE. If BITPIX is
+	# not one of the MII types then rft_read_pixels returns an ERROR.
+
+	call amovkl (long(1), v, IM_MAXDIM)
+	switch (BITPIX(fits)) {
+	case FITS_REAL:
+
+	    # Allocate temporary space.
+	    call malloc (tempbuf, npix, TY_REAL)
+
+	    # Initialize the read.
+	    i = rft_init_read_pixels (npix_record, BITPIX(fits), LSBF, TY_REAL)
+
+	    # Turn on the ieee NaN mapping.
+	    call ieesnanr (blank)
+	    #call ieemapr (YES, NO)
+	    #call ieezstatr ()
+	    NBPIX(im) = 0
+
+	    # Allocate the space for the output line, read in the image
+	    # line, convert from the ieee to native format, and compute the
+	    # minimum and maximum.
+
+	    do il = 1, nlines {
+	        call rft_put_image_line (im, buf, v, PIXTYPE(im))
+	        if (rft_read_pixels (fits_fd, Memr[tempbuf], npix,
+	            NRECORDS(fits), blksize) != npix)
+		    call printf ("Error reading FITS data\n")
+	        if (SCALE(fits) == YES)
+	            call rft_rscale_pix (Memr[tempbuf], buf, npix,
+		        FITS_BSCALE(fits), FITS_BZERO(fits), PIXTYPE(im))
+		else
+		    call rft_rchange_pix (Memr[tempbuf], buf, npix, PIXTYPE(im))
+	        call rft_pix_limits (buf, npix, PIXTYPE(im), linemin, linemax)
+	        lirafmax = max (lirafmax, linemax)
+	        lirafmin = min (lirafmin, linemin)
+	    }
+
+	    # Set the number of bad pixels.
+	    call ieestatr (NBPIX(im), ndummy)
+
+	    # Free space.
+	    call mfree (tempbuf, TY_REAL)
+
+	case FITS_DOUBLE:
+
+	    # Allocate temporary space.
+	    call malloc (tempbuf, npix, TY_DOUBLE)
+
+	    # Initialize the read.
+	    i = rft_init_read_pixels (npix_record, BITPIX(fits), LSBF,
+	        TY_DOUBLE)
+
+	    # Turn on the ieee NaN mapping.
+	    dblank = blank
+	    call ieesnand (dblank)
+	    #call ieemapd (YES, NO)
+	    #call ieezstatd ()
+	    NBPIX(im) = 0
+
+	    # Allocate the space for the output line, read in the image
+	    # line, convert from the ieee to native format, and compute the
+	    # minimum and maximum.
+
+	    do il = 1, nlines {
+	        call rft_put_image_line (im, buf, v, PIXTYPE(im))
+	        if (rft_read_pixels (fits_fd, Memd[tempbuf], npix,
+	            NRECORDS(fits), blksize) != npix)
+		    call printf ("Error reading FITS data\n")
+	        if (SCALE(fits) == YES)
+	            call rft_dscale_pix (Memd[tempbuf], buf, npix,
+		        FITS_BSCALE(fits), FITS_BZERO(fits), PIXTYPE(im))
+		else
+		    call rft_dchange_pix (Memd[tempbuf], buf, npix, PIXTYPE(im))
+	        call rft_pix_limits (buf, npix, PIXTYPE(im), linemin, linemax)
+		if (IS_INDEFR(linemax))
+		    lirafmax = INDEFR
+		else
+	            lirafmax = max (lirafmax, linemax)
+		if (IS_INDEFR(linemin))
+		    lirafmin = INDEFR
+		else
+	            lirafmin = min (lirafmin, linemin)
+	    }
+
+	    # Set the number of bad pixels.
+	    call ieestatd (NBPIX(im), ndummy)
+
+	    # Free space.
+	    call mfree (tempbuf, TY_DOUBLE)
+
+	default:
+
+	    # Allocate the required space.
+	    call malloc (tempbuf, npix, TY_LONG)
+
+	    # Allocate the space for the output line, read in the image
+	    # line, convert from the ieee to native format, and compute the
+	    # minimum and maximum.
+
+	    i = rft_init_read_pixels (npix_record, BITPIX(fits), LSBF, TY_LONG)
+	    do il = 1, nlines {
+	        call rft_put_image_line (im, buf, v, PIXTYPE(im))
+	        if (rft_read_pixels (fits_fd, Meml[tempbuf], npix,
+	            NRECORDS(fits), blksize) != npix)
+		    call printf ("Error reading FITS data\n")
+	        if (SCALE(fits) == YES)
+	            call rft_lscale_pix (Meml[tempbuf], buf, npix,
+		        FITS_BSCALE(fits), FITS_BZERO(fits), PIXTYPE(im))
+	        else
+		    call rft_lchange_pix (Meml[tempbuf], buf, npix, PIXTYPE(im))
+	        if (BLANKS(fits) == YES)
+	            call rft_map_blanks (Meml[tempbuf], buf, npix, PIXTYPE(im),
+		        BLANK_VALUE(fits), blank, NBPIX(im))
+	        call rft_pix_limits (buf, npix, PIXTYPE(im), linemin, linemax)
+	        lirafmax = max (lirafmax, linemax)
+	        lirafmin = min (lirafmin, linemin)
+	    }
+
+	    # Free space.
+	    call mfree (tempbuf, TY_LONG)
+	}
+
+	IRAFMIN(im) = lirafmin
+	IRAFMAX(im) = lirafmax
+	LIMTIME(im) = clktime (long(0))
+
+	if (NBPIX (im) != 0) {
+	    call printf ("Warning: %d bad pixels replaced in image\n")
+		call pargl (NBPIX (im))
+	}
+	if (IS_INDEFR(lirafmax) || lirafmax > MAX_REAL) {
+	    call printf ("Warning: image contains pixel values > %g\n")
+		call pargr (MAX_REAL)
+	}
+	if (IS_INDEFR(lirafmin) || lirafmin < -MAX_REAL) {
+	    call printf ("Warning: image contains pixel values < %g\n")
+		call pargr (-MAX_REAL)
+	}
+end
+
+
+# RFT_SET_IMAGE_HEADER -- Set remaining header fields not set in
+# rft_read_header.
+
+procedure rft_set_image_header (fits, im)
+
+pointer	fits		# FITS data structure
+pointer	im		# IRAF image pointer
+
+include	"rfits.com"
+
+begin
+	# Determine data type from BITPIX if user data type not specified.
+
+	if (data_type == ERR) {
+	    if (BITPIX(fits) < 0) {
+		if (abs (BITPIX(fits)) <= (SZ_REAL * SZB_CHAR * NBITS_BYTE))
+		    PIXTYPE(im) = TY_REAL
+		else
+		    PIXTYPE(im) = TY_DOUBLE
+	    } else if (SCALE(fits) == YES) {
+		PIXTYPE(im) = TY_REAL
+	    } else {
+	        if (BITPIX(fits) <= (SZ_SHORT * SZB_CHAR * NBITS_BYTE))
+		    PIXTYPE(im) = TY_SHORT
+		else
+		    PIXTYPE(im) = TY_LONG
+	    }
+
+	} else
+	    PIXTYPE(im) = data_type
+end
+
+
+# RFT_SET_PRECISION -- Procedure to determine the precision of the FITS data
+# type.
+
+procedure rft_set_precision (bitpix, precision)
+
+int	bitpix			# FITS bits per pixel
+int	precision		# FITS decimal digits of precision
+
+begin
+	switch (bitpix) {
+	case FITS_BYTE:
+	    precision = FITSB_PREC
+	case FITS_SHORT:
+	    precision = FITSS_PREC
+	case FITS_LONG:
+	    precision = FITSL_PREC
+	default:
+	    call error (16, "RFT_SET_PRECISION: Unknown FITS type")
+	}
+end
+
+
+# RFT_PUT_IMAGE_LINE -- Procedure to output an image line to and IRAF file.
+
+procedure rft_put_image_line (im, buf, v, data_type)
+
+pointer	im			# IRAF image descriptor
+pointer	buf			# Pointer to output image line
+long	v[ARB]			# imio pointer
+int	data_type		# output pixel type
+
+int	impnll(), impnlr(), impnld(), impnlx()
+errchk	impnll, impnlr, impnld, impnlx
+
+begin
+	switch (data_type) {
+	case TY_SHORT, TY_INT, TY_USHORT, TY_LONG:
+	    if (impnll (im, buf, v) == EOF)
+		call error (3, "RFT_PUT_IMAGE_LINE: Error writing FITS data")
+	case TY_REAL:
+	    if (impnlr (im, buf, v) == EOF)
+		call error (3, "RFT_PUT_IMAGE_LINE: Error writing FITS data")
+	case TY_DOUBLE:
+	    if (impnld (im, buf, v) == EOF)
+		call error (3, "RFT_PUT_IMAGE_LINE: Error writing FITS data")
+	case TY_COMPLEX:
+	    if (impnlx (im, buf, v) == EOF)
+		call error (3, "RFT_PUT_IMAGE_LINE: Error writing FITS data")
+	default:
+	    call error (10, "RFT_PUT_IMAGE_LINE: Unsupported IRAF image type")
+	}
+end
+
+
+# RFT_RSCALE_PIX -- Procedure to convert an IRAF image line from type real
+# to the requested output data type with optional scaling using the
+# FITS parameters BSCALE and BZERO.
+
+procedure rft_rscale_pix (inbuf, outbuf, npix, bscale, bzero, data_type)
+
+real	inbuf[ARB]		# buffer of FITS integers
+pointer	outbuf			# pointer to output image line
+int	npix			# number of pixels
+double	bscale, bzero		# FITS bscale and bzero
+int	data_type		# IRAF image pixel type
+
+errchk	altmdr, achtrl, amovr, achtrd, achtrx
+
+begin
+	switch (data_type) {
+	case TY_SHORT, TY_USHORT, TY_INT, TY_LONG:
+	    call altmdr (inbuf, inbuf, npix, bscale, bzero)
+	    call achtrl (inbuf, Meml[outbuf], npix)
+	case TY_REAL:
+	    call altmdr (inbuf, inbuf, npix, bscale, bzero)
+	    call amovr (inbuf, Memr[outbuf], npix)
+	case TY_DOUBLE:
+	    call altmdr (inbuf, inbuf, npix, bscale, bzero)
+	    call achtrd (inbuf, Memd[outbuf], npix)
+	case TY_COMPLEX:
+	    call altmdr (inbuf, inbuf, npix, bscale, bzero)
+	    call achtrx (inbuf, Memx[outbuf], npix)
+	default:
+	    call error (10, "RFT_SCALE_LINE: Illegal IRAF image type")
+	}
+end
+
+
+# RFT_DSCALE_PIX -- Procedure to convert an IRAF image line from type double
+# to the requested output data type with optional scaling using the
+# FITS parameters BSCALE and BZERO.
+
+procedure rft_dscale_pix (inbuf, outbuf, npix, bscale, bzero, data_type)
+
+double	inbuf[ARB]		# buffer of FITS integers
+pointer	outbuf			# pointer to output image line
+int	npix			# number of pixels
+double	bscale, bzero		# FITS bscale and bzero
+int	data_type		# IRAF image pixel type
+
+errchk	altmd, achtdl, amovd, achtdr, achtdx
+
+begin
+	switch (data_type) {
+	case TY_SHORT, TY_USHORT, TY_INT, TY_LONG:
+	    call altmd (inbuf, inbuf, npix, bscale, bzero)
+	    call achtdl (inbuf, Meml[outbuf], npix)
+	case TY_REAL:
+	    call altmd (inbuf, inbuf, npix, bscale, bzero)
+	    call achtdr (inbuf, Memr[outbuf], npix)
+	case TY_DOUBLE:
+	    call altmd (inbuf, inbuf, npix, bscale, bzero)
+	    call amovd (inbuf, Memd[outbuf], npix)
+	case TY_COMPLEX:
+	    call altmd (inbuf, inbuf, npix, bscale, bzero)
+	    call achtdx (inbuf, Memx[outbuf], npix)
+	default:
+	    call error (10, "RFT_SCALE_LINE: Illegal IRAF image type")
+	}
+end
+
+
+
+# RFT_LSCALE_PIX -- Procedure to convert an IRAF image line from type long
+# to the requested output data type with optional scaling using the
+# FITS parameters BSCALE and BZERO.
+
+procedure rft_lscale_pix (inbuf, outbuf, npix, bscale, bzero, data_type)
+
+long	inbuf[ARB]		# buffer of FITS integers
+pointer	outbuf			# pointer to output image line
+int	npix			# number of pixels
+double	bscale, bzero		# FITS bscale and bzero
+int	data_type		# IRAF image pixel type
+
+errchk	achtll, achtlr, achtld, achtlx
+errchk	altml, altmr, altmd, altmx
+
+begin
+	switch (data_type) {
+	case TY_SHORT, TY_USHORT, TY_INT, TY_LONG:
+	    call achtll (inbuf, Meml[outbuf], npix)
+	    call altml (Meml[outbuf], Meml[outbuf], npix, bscale, bzero)
+	case TY_REAL:
+	    call altmlr (inbuf, Memr[outbuf], npix, bscale, bzero)
+	    #call achtlr (inbuf, Memr[outbuf], npix)
+	    #call altmdr (Memr[outbuf], Memr[outbuf], npix, bscale, bzero)
+	case TY_DOUBLE:
+	    call achtld (inbuf, Memd[outbuf], npix)
+	    call altmd (Memd[outbuf], Memd[outbuf], npix, bscale, bzero)
+	case TY_COMPLEX:
+	    call achtlx (inbuf, Memx[outbuf], npix)
+	    call altmx (Memx[outbuf], Memx[outbuf], npix, real (bscale),
+	        real (bzero))
+	default:
+	    call error (10, "RFT_SCALE_LINE: Illegal IRAF image type")
+	}
+end
+
+
+# RFT_RCHANGE_PIX -- Procedure to change a line of real numbers to the
+# IRAF image type.
+
+procedure rft_rchange_pix (inbuf, outbuf, npix, data_type)
+
+real	inbuf[ARB]		# array of FITS integers
+pointer	outbuf			# pointer to IRAF image line
+int	npix			# number of pixels
+int	data_type		# IRAF pixel type
+
+errchk	achtrl, amovr, achtrd, achtrx
+
+begin
+	switch (data_type) {
+	case TY_SHORT, TY_INT, TY_USHORT, TY_LONG:
+	    call achtrl (inbuf, Meml[outbuf], npix)
+	case TY_REAL:
+	    call amovr (inbuf, Memr[outbuf], npix)
+	case TY_DOUBLE:
+	    call achtrd (inbuf, Memd[outbuf], npix)
+	case TY_COMPLEX:
+	    call achtrx (inbuf, Memx[outbuf], npix)
+	default:
+	    call error (10, "RFT_RCHANGE_LINE: Illegal IRAF image type")
+	}
+end
+
+
+# RFT_DCHANGE_PIX -- Procedure to change a line of double precision numbers
+# to the IRAF image type.
+
+procedure rft_dchange_pix (inbuf, outbuf, npix, data_type)
+
+double  inbuf[ARB]		# array of FITS integers
+pointer	outbuf			# pointer to IRAF image line
+int	npix			# number of pixels
+int	data_type		# IRAF pixel type
+
+errchk	achtdl, achtdr, amovd, achtdx
+
+begin
+	switch (data_type) {
+	case TY_SHORT, TY_INT, TY_USHORT, TY_LONG:
+	    call achtdl (inbuf, Meml[outbuf], npix)
+	case TY_REAL:
+	    call achtdr (inbuf, Memr[outbuf], npix)
+	case TY_DOUBLE:
+	    call amovd (inbuf, Memd[outbuf], npix)
+	case TY_COMPLEX:
+	    call achtdx (inbuf, Memx[outbuf], npix)
+	default:
+	    call error (10, "RFT_DCHANGE_LINE: Illegal IRAF image type")
+	}
+end
+
+
+
+# RFT_LCHANGE_PIX -- Procedure to change a line of long integers to the
+# IRAF image type.
+
+procedure rft_lchange_pix (inbuf, outbuf, npix, data_type)
+
+long	inbuf[ARB]		# array of FITS integers
+pointer	outbuf			# pointer to IRAF image line
+int	npix			# number of pixels
+int	data_type		# IRAF pixel type
+
+begin
+	switch (data_type) {
+	case TY_SHORT, TY_USHORT, TY_INT, TY_LONG:
+	    call achtll (inbuf, Meml[outbuf], npix)
+	case TY_REAL:
+	    call achtlr (inbuf, Memr[outbuf], npix)
+	case TY_DOUBLE:
+	    call achtld (inbuf, Memd[outbuf], npix)
+	case TY_COMPLEX:
+	    call achtlx (inbuf, Memx[outbuf], npix)
+	default:
+	    call error (10, "RFT_CHANGE_LINE: Illegal IRAF image type")
+	}
+end
+
+
+# RFT_MAP_BLANKS -- Map the blank pixels.  Currently only the number of blank
+# pixels is determined without an further mapping.
+
+procedure rft_map_blanks (a, buf, npts, pixtype, blank_value, blank, nbadpix)
+
+long	a[ARB]		# integer input buffer
+pointer	buf		# pointer to output image buffer
+int	npts		# number of points
+int	pixtype		# image data type
+long	blank_value	# FITS blank value
+real	blank		# user blank value
+long	nbadpix		# number of bad pixels
+
+int	i
+
+begin
+	# Do blank mapping here
+	switch (pixtype) {
+	case TY_SHORT, TY_INT, TY_USHORT, TY_LONG:
+	    do i = 1, npts {
+	        if (a[i] == blank_value) {
+		    nbadpix = nbadpix + 1
+		    Meml[buf+i-1] = blank
+		}
+	    }
+	case TY_REAL:
+	    do i = 1, npts {
+	        if (a[i] == blank_value) {
+		    nbadpix = nbadpix + 1
+		    Memr[buf+i-1] = blank
+		}
+	    }
+	case TY_DOUBLE:
+	    do i = 1, npts {
+	        if (a[i] == blank_value) {
+		    nbadpix = nbadpix + 1
+		    Memd[buf+i-1] = blank
+		}
+	    }
+	case TY_COMPLEX:
+	    do i = 1, npts {
+	        if (a[i] == blank_value) {
+		    nbadpix = nbadpix + 1
+		    Memx[buf+i-1] = blank
+		}
+	    }
+	}
+end
+
+
+# RFT_PIX_LIMITS -- Procedure to determine to maxmimum and minimum values in a
+# line. Note that double precision is somewhat of a special case because
+# MAX_DOUBLE is a lot less than the maximum permitted ieee numbers for iraf.
+
+procedure rft_pix_limits (buf, npix, pixtype, linemin, linemax)
+
+pointer	buf				# pointer to IRAF image line
+int	npix				# number of pixels
+int	pixtype				# output data type
+real	linemax, linemin		# min and max pixel values
+
+long	lmax, lmin
+real	rmax, rmin
+double	dmax, dmin
+complex	xmax, xmin
+
+begin
+	switch (pixtype) {
+	case TY_SHORT, TY_INT, TY_USHORT, TY_LONG:
+	    call aliml (Meml[buf], npix, lmin, lmax)
+	    linemax = lmax
+	    linemin = lmin
+	case TY_REAL:
+	    call alimr (Memr[buf], npix, rmin, rmax)
+	    linemax = rmax
+	    linemin = rmin
+	case TY_DOUBLE:
+	    call alimd (Memd[buf], npix, dmin, dmax)
+	    if (dmax > MAX_REAL)
+	        linemax = INDEFR
+	    else
+		linemax = dmax
+	    if (dmin < -MAX_REAL)
+		linemin = INDEFR
+	    else
+	        linemin = dmin
+	case TY_COMPLEX:
+	    call alimx (Memx[buf], npix, xmin, xmax)
+	    linemax = xmax
+	    linemin = xmin
+	default:
+	    call error (30, "RFT_PIX_LIMITS: Unknown IRAF type")
+	}
+end
+
+
+# ALTMDR -- procedure to scale a long vector into a real vector using
+# double precision constants to preserve accuracy
+
+procedure altmlr (a, b, npix, bscale, bzero)
+
+long	a[ARB]		# input array
+real	b[ARB]		# output array
+int	npix		# number of pixels
+double	bscale, bzero	# scaling parameters
+
+int	i
+
+begin
+	do i = 1, npix
+	    b[i] = a[i] * bscale + bzero
+end
+
+
+# ALTMDR -- procedure to scale a real vector with double precision constants.
+
+procedure altmdr (a, b, npix, bscale, bzero)
+
+real	a[ARB]		# input array
+real	b[ARB]		# output array
+int	npix		# number of pixels
+double	bscale, bzero	# scaling parameters
+
+int	i
+
+begin
+	do i = 1, npix
+	    b[i] = a[i] * bscale + bzero
+end
diff --git a/pkg/obsolete/fits/fits_rpixels.x b/pkg/obsolete/fits/fits_rpixels.x
new file mode 100644
index 00000000..dfce821d
--- /dev/null
+++ b/pkg/obsolete/fits/fits_rpixels.x
@@ -0,0 +1,154 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include <fset.h>
+include <mii.h>
+include <mach.h>
+
+# RFT_INIT_READ_PIXELS and READ_PIXELS -- Read pixel data with record buffering
+# and data type conversion.  The input data must meet the MII standard
+# except for possibly in the case of integers having the least significant
+# byte first.
+#
+# Read data in records of len_record and convert to the specified IRAF
+# data type.  Successive calls of rft_read_pixels returns the next npix pixels.
+# Read_pixels returns EOF or the number of pixels converted.
+# Init_read_pixels must be called before read_pixels.
+#
+# Error conditions are:
+# 1. A short input record
+# 2. Error in converting the pixels by miiup.
+#
+# This routine is based on the MII unpack routine which is machine dependent.
+# The bitpix must correspond to an MII type.  If the lsbf (least significant
+# byte first) flag is YES then the pixels do not satisfy the MII standard.
+# In this case the bytes are first swapped into most significant byte first
+# before the MII unpack routine is called.
+
+int procedure rft_init_read_pixels (npix_record, bitpix, lsbf, spp_type)
+
+int	npix_record	# number of pixels per input record
+int	bitpix		# bits per pixel (must correspond to an MII type)
+int	lsbf		# byte swap?
+int	spp_type	# SPP data type to be returned
+
+# entry rft_read_pixels (fd, buffer, npix)
+
+int	rft_read_pixels
+int	fd		# input file descriptor
+char	buffer[1]	# output buffer
+int	npix		# number of pixels to read
+
+int	swap
+int	ty_mii, ty_spp, npix_rec, nch_rec, sz_rec, nchars, len_mii, recptr
+int	bufsize, i, n, ip, op
+pointer	mii, spp
+
+int	rft_getbuf(), sizeof(), miilen()
+errchk	miilen, mfree, malloc, rft_getbuf, miiupk
+data	mii/NULL/, spp/NULL/
+
+begin
+	ty_mii = bitpix
+	ty_spp = spp_type
+	swap = lsbf
+	npix_rec = npix_record
+	nch_rec = npix_rec * sizeof (ty_spp)
+
+	len_mii = miilen (npix_rec, ty_mii)
+	sz_rec = len_mii * SZ_INT
+
+	if (mii != NULL)
+	    call mfree (mii, TY_INT)
+	call malloc (mii, len_mii, TY_INT)
+
+	if (spp != NULL)
+	    call mfree (spp, TY_CHAR)
+	call malloc (spp, nch_rec, TY_CHAR)
+
+	ip = nch_rec
+	return (OK)
+
+entry	rft_read_pixels (fd, buffer, npix, recptr, bufsize)
+
+	nchars = npix * sizeof (ty_spp)
+	op = 0
+
+	repeat {
+
+	    # If data is exhausted read the next record
+	    if (ip == nch_rec) {
+
+		i = rft_getbuf (fd, Memi[mii], sz_rec, bufsize, recptr)
+		if (i == EOF)
+		    return (EOF)
+
+		if (swap == YES)
+		    switch (ty_mii) {
+		    case MII_SHORT:
+			call bswap2 (Memi[mii], 1, Memi[mii], 1,
+				sz_rec * SZB_CHAR)
+		    case MII_LONG:
+			call bswap4 (Memi[mii], 1, Memi[mii], 1,
+				sz_rec * SZB_CHAR)
+		    }
+
+		call miiupk (Memi[mii], Memc[spp], npix_rec, ty_mii, ty_spp)
+
+		ip = 0
+		#recptr = recptr + 1
+	    }
+
+	    n = min (nch_rec - ip, nchars - op)
+	    call amovc (Memc[spp+ip], buffer[1+op], n)
+	    ip = ip + n
+	    op = op + n
+
+	} until (op == nchars)
+
+	return (npix)
+end
+
+
+# RFT_GETBUF -- Procedure to get the buffer.
+
+int procedure rft_getbuf (fd, buf, sz_rec, bufsize, recptr)
+
+int	fd		# file descriptor
+char	buf[ARB]	# buffer to be filled
+int	sz_rec		# size in chars of record to be read
+int	bufsize		# buffer size in records
+int	recptr		# last successful FITS record read
+
+int	i, nchars
+int	read(), fstati()
+errchk	read
+
+begin
+	nchars = 0
+	repeat {
+	    iferr {
+	        i = read (fd, buf[nchars+1], sz_rec - nchars)
+	    } then {
+	        call printf ("Error reading FITS record %d\n")
+	        if (mod (recptr + 1, bufsize) == 0)
+		    call pargi ((recptr + 1) / bufsize)
+	        else
+		    call pargi ((recptr + 1) / bufsize + 1)
+	        call fseti (fd, F_VALIDATE, fstati (fd, F_SZBBLK) / SZB_CHAR)
+	        i = read (fd, buf[nchars+1], sz_rec - nchars)
+	    }
+
+	    if (i == EOF)
+		break
+	    else
+	        nchars = nchars + i
+
+	} until (nchars >= sz_rec)
+
+	if ((i == EOF) && (nchars == 0))
+	    return (EOF)
+	else {
+	    recptr = recptr + 1
+	    return (nchars)
+	}
+end
diff --git a/pkg/obsolete/fits/fits_wheader.x b/pkg/obsolete/fits/fits_wheader.x
new file mode 100644
index 00000000..f00a0de6
--- /dev/null
+++ b/pkg/obsolete/fits/fits_wheader.x
@@ -0,0 +1,471 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include <imhdr.h>
+include <mach.h>
+include <fset.h>
+include "wfits.h"
+
+# WFT_WRITE_HEADER -- Write the FITS headers. The FITS header
+# parameters are encoded one by one until the FITS END keyword is detected.
+# If the long_header switch is set the full FITS header is printed on the
+# standard output. If the short header parameter is specified only the image
+# title and dimensions are printed.
+
+procedure wft_write_header (im, fits, fits_fd)
+
+pointer	im		# pointer to the IRAF image
+pointer	fits		# pointer to the FITS structure
+int	fits_fd		# the FITS file descriptor
+
+char	card[LEN_CARD+1], trim_card[LEN_CARD+1]
+int	nrecords, recntr, cardptr, cardcnt, stat, cards_per_rec, i
+int	wft_card_encode(), wft_set_bitpix(), sizeof(), strncmp()
+int	wft_init_card_encode(), fstati()
+
+errchk	wft_set_bitpix, wft_get_iraf_typestring, wft_set_scale, wft_set_blank
+errchk	wft_fits_set_scale, wft_init_card_encode, wft_card_encode
+errchk	wft_init_write_pixels, wft_write_pixels, wft_write_last_record
+
+include "wfits.com"
+
+begin
+	# SET the data and FITS bits per pixel.
+
+	DATA_BITPIX(fits) = sizeof (PIXTYPE(im)) * SZB_CHAR * NBITS_BYTE
+	FITS_BITPIX(fits) = wft_set_bitpix (bitpix, PIXTYPE(im),
+	    DATA_BITPIX(fits))
+
+	# Calculate the FITS bscale and bzero parameters. Notice for the
+	# time being that scaling is turned off if IEEE floating point
+	# output is selected. May decide to change this later after
+	# checking the specifications.
+
+	if (FITS_BITPIX(fits) < 0) {
+
+	    IRAFMIN(fits) = IM_MIN(im)
+	    IRAFMAX(fits) = IM_MAX(im)
+	    SCALE(fits) = NO
+	    BZERO(fits) = 0.0d0
+	    BSCALE(fits) = 1.0d0
+
+	} else if (autoscale == YES) {
+
+	    call wft_get_tape_limits (FITS_BITPIX(fits), TAPEMIN(fits),
+	        TAPEMAX(fits))
+	    call wft_data_limits (im, IRAFMIN(fits), IRAFMAX(fits))
+	    call wft_fits_set_scale (im, DATA_BITPIX(fits), FITS_BITPIX(fits),
+	        IRAFMIN(fits), IRAFMAX(fits), TAPEMIN(fits), TAPEMAX(fits),
+		SCALE(fits), BSCALE(fits), BZERO(fits))
+
+	} else {
+
+	    IRAFMIN(fits) = IM_MIN(im)
+	    IRAFMAX(fits) = IM_MAX(im)
+	    SCALE(fits) = scale
+	    BZERO(fits) = bzero
+	    BSCALE(fits) = bscale
+	}
+
+	# If blanks in the image set the blank parameter. Currently information
+	# on blanks is not written out so this is effectively a null operation
+	# in IRAF.
+
+	if (NBPIX(im) > 0)
+	    call wft_set_blank (FITS_BITPIX(fits), BLANK(fits),
+	        BLANK_STRING(fits))
+
+	# Set the IRAF datatype parameter.
+	call wft_get_iraf_typestring (PIXTYPE(im), TYPE_STRING(fits))
+
+	# Initialize the card counters. These counters are used only for
+	# information printed to the standard output.
+
+	recntr = 1
+	cardptr = 1
+	cardcnt = 1
+	cards_per_rec = len_record / LEN_CARD
+
+	# Get set up to write the FITS header. Initialize for an ASCII write.
+	stat = wft_init_card_encode (im, fits)
+	if (make_image == YES)
+	    call wft_init_wrt_pixels (len_record, TY_CHAR, FITS_BYTE, blkfac)
+
+	# Print short header.
+	if (short_header == YES && long_header == NO) {
+
+	    call printf ("%-20.20s  ")
+		call pargstr (OBJECT(im))
+	    do i = 1, NAXIS(im) {
+		if (i == 1) {
+		    call printf ("Size = %d")
+			call pargl (NAXISN(im,i))
+		} else {
+		    call printf (" x %d")
+			call pargl (NAXISN(im,i))
+		}
+	    }
+	    call printf ("\n")
+
+	    call strlwr (TYPE_STRING(fits))
+	    call printf ("\tpixtype=%s bitpix=%d")
+		call pargstr (TYPE_STRING(fits))
+		call pargi (FITS_BITPIX(fits))
+
+	    if (fstati (fits_fd, F_BLKSIZE) == 0) {
+		call printf (" blkfac=%d")
+		    call pargi (blkfac)
+	    } else
+		call printf (" blkfac=fixed")
+
+	    if (SCALE(fits) == YES) {
+	        call printf (" bscale=%.7g bzero=%.7g\n")
+		    call pargd (BSCALE(fits))
+		    call pargd (BZERO(fits))
+	    } else
+		call printf (" scaling=none\n")
+	    call strupr (TYPE_STRING(fits))
+	}
+
+	# Write the cards to the FITS header.
+	repeat {
+
+	    # Encode the card.
+	    stat = wft_card_encode (im, fits, card)
+	    if (stat == NO)
+		next
+
+	    # Write the card to the output file if make_image is yes.
+	    if (make_image == YES)
+	        call wft_write_pixels (fits_fd, card, LEN_CARD)
+
+	    # Trim the card and write is to the standard output if 
+	    # long_header is yes.
+
+	    if (long_header == YES) {
+		call wft_trimstr (card, trim_card, LEN_CARD)
+	        call printf ("%2d/%2d:--  %s\n")
+		    call pargi (recntr)
+		    call pargi (cardptr)
+		    call pargstr (trim_card)
+	    }
+
+	    if (mod (cardcnt, cards_per_rec) == 0) {
+	        recntr = recntr + 1
+	        cardptr = 1
+	    } else
+		cardptr = cardptr + 1
+	    cardcnt = cardcnt + 1
+
+	} until (strncmp (card, "END     ", LEN_KEYWORD) == 0)
+
+	# Issue warning about possible precision loss. Comment this out
+	# for time being, since the short header was modified.
+	#if (SCALE(fits) == YES && bitpix != ERR) {
+	    #call printf (
+	    #"\tDefault bitpix overridden: maximum precision loss ~%.7g\n")
+		#call pargd (BSCALE(fits))
+	#}
+
+	# Write the last header records.
+	if (make_image == YES) {
+	    call wft_write_last_record (fits_fd, nrecords)
+	    if (short_header == YES || long_header == YES) {
+	        call printf ("\t%d Header  ")
+		    call pargi (nrecords)
+	    }
+	}
+end
+
+
+# WFT_SET_BITPIX -- This procedure sets the FITS bitpix for each image based on
+# either the user given value or the precision of the IRAF data. Notice that
+# the user must explicitly set the bitpix parameter to -16 or -32 to select
+# the IEEE output format.
+
+int procedure wft_set_bitpix (bitpix, datatype, data_bitpix)
+
+int	bitpix		# the user set bits per pixel, ERR or legal bitpix
+int	datatype	# the IRAF image data type
+int	data_bitpix	# the bits per pixel in the data
+
+begin
+	if (bitpix == ERR) {
+	    switch (datatype) {
+	    case TY_SHORT, TY_INT, TY_USHORT, TY_LONG:
+	        if (data_bitpix <= FITS_BYTE)
+		    return (FITS_BYTE)
+	        else if (data_bitpix <= FITS_SHORT) {
+		    if (datatype == TY_USHORT)
+			return (FITS_LONG)
+		    else
+		        return (FITS_SHORT)
+	        } else
+		    return (FITS_LONG)
+	    case TY_REAL, TY_COMPLEX:
+		return (FITS_REAL)
+	    case TY_DOUBLE:
+		 return (FITS_DOUBLE)
+	    default:
+		call error (2, "SET_BITPIX: Unknown IRAF data type.")
+	    }
+	} else
+	    return (bitpix)
+end
+
+
+# WFT_GET_IRAF_TYPESTRING -- Procedure to set the iraf datatype keyword.
+# Permitted strings are INTEGER, FLOATING or COMPLEX.
+
+procedure wft_get_iraf_typestring (datatype, type_str)
+
+int	datatype	# the IRAF data type
+char	type_str[ARB]	# the output IRAF type string
+
+begin
+	switch (datatype) {
+	case TY_SHORT:
+	    call strcpy ("SHORT", type_str, LEN_STRING) 
+	case TY_USHORT:
+	    call strcpy ("USHORT", type_str, LEN_STRING) 
+	case TY_INT: 
+	    call strcpy ("INTEGER", type_str, LEN_STRING) 
+	case TY_LONG: 
+	    call strcpy ("LONG", type_str, LEN_STRING) 
+	case TY_REAL:
+	    call strcpy ("REAL", type_str, LEN_STRING)
+	case TY_DOUBLE:
+	    call strcpy ("DOUBLE", type_str, LEN_STRING)
+	case TY_COMPLEX:
+	    call strcpy ("COMPLEX", type_str, LEN_STRING)
+	default:
+	    call error (3, "IRAF_TYPE: Unknown IRAF image type.")
+	}
+end
+
+
+# WFT_FITS_SET_SCALE -- Procedure to set the FITS scaling parameters if
+# autoscaling is enabled.
+
+procedure wft_fits_set_scale (im, data_bitpix, fits_bitpix, irafmin, irafmax,
+	tapemin, tapemax, scale, bscale, bzero )
+
+pointer	im		# pointer to IRAF image
+int	data_bitpix	# bits per pixel of data
+int	fits_bitpix	# fits bits per pixel
+real	irafmin		# minimum picture value
+real	irafmax		# maximum picture value
+double	tapemin		# minimum tape value
+double	tapemax		# maximum tape value
+int	scale		# scale data ?
+double	bscale		# FITS bscale
+double	bzero		# FITS bzero
+
+errchk	wft_set_scale
+
+begin
+	switch (PIXTYPE(im)) {
+	case TY_SHORT, TY_INT, TY_LONG:
+	    if (data_bitpix > fits_bitpix) {
+		scale = YES
+		call wft_set_scale (fits_bitpix, irafmin, irafmax, tapemin, 
+		    tapemax, bscale, bzero)
+	    } else {
+		scale = NO
+		bscale = 1.0d0
+		bzero = 0.0d0
+	    }
+	case TY_USHORT:
+	    if (data_bitpix > fits_bitpix) {
+		scale = YES
+		call wft_set_scale (fits_bitpix, irafmin, irafmax, tapemin, 
+		    tapemax, bscale, bzero)
+	    } else if (data_bitpix == fits_bitpix) {
+		scale = YES
+		bscale = 1.0d0
+		bzero = 3.2768d4
+	    } else {
+		scale = NO
+		bscale = 1.0d0
+		bzero = 0.0d0
+	    }
+	case TY_REAL, TY_DOUBLE, TY_COMPLEX:
+	    scale = YES
+	    call wft_set_scale (fits_bitpix, irafmin, irafmax, tapemin, tapemax,
+	        bscale, bzero)
+	default:
+	    call error (1, "WRT_HEADER: Unknown IRAF image type.")
+	}
+
+end
+
+
+# WFT_SET_SCALE -- This procedure calculates bscale and bzero for each frame
+# from the FITS bitpix and the maximum and minimum data values.
+
+procedure wft_set_scale (fits_bitpix, datamin, datamax, mintape, maxtape,
+	bscale, bzero)
+
+int	fits_bitpix		# the FITS integer bits per pixels
+real	datamax, datamin	# the IRAF image data minimum and maximum
+double	mintape, maxtape	# min and max FITS tape values
+double	bscale, bzero		# the calculated bscale and bzero values
+
+double	maxdata, mindata, num, denom
+bool	rft_equald()
+
+begin
+	# Calculate the maximum and minimum values in the data.
+	maxdata = datamax + abs ((datamax / (10.0 ** (NDIGITS_RP - 1))))
+	mindata = datamin - abs ((datamin / (10.0 ** (NDIGITS_RP - 1))))
+	denom =  maxtape - mintape
+	num = maxdata - mindata
+	#denom = denom - denom / (1.0d1 ** (NDIGITS_RP - 1))
+	#num = num + num / (1.0d1 ** (NDIGITS_RP - 1))
+
+	# Check for constant image case.
+	mindata = datamin
+	maxdata = datamax
+	if (rft_equald (num, 0.0d0)) {
+	    bscale = 1.0d0
+	    bzero = maxdata
+	} else {
+	    bscale = num / denom
+	    bzero = (maxtape / denom) * mindata - (mintape / denom) * maxdata
+	}
+end
+
+
+# WFT_GET_TAPE_LIMITS -- Procedure for calculating the maximum and minimum FITS
+# integer values from the FITS bitpix.
+
+procedure wft_get_tape_limits (fits_bitpix, mintape, maxtape)
+
+int	fits_bitpix		# the bits per pixel of a FITS integer
+double	maxtape, mintape	# the maximun and minimum FITS tape integers
+
+begin
+	switch (fits_bitpix) {
+	case FITS_BYTE:
+	    maxtape = BYTE_MAX
+	    mintape = BYTE_MIN
+	case FITS_SHORT:
+	    maxtape = SHORT_MAX
+	    mintape = SHORT_MIN
+	case FITS_LONG:
+	    maxtape = LONG_MAX
+	    mintape = LONG_MIN
+	default:
+	    call error (4, "TAPE_LIMITS: Unknown FITS type.")
+	}
+end
+
+
+# WFT_SET_BLANK -- Determine the FITS integer value for a blank pixel from the
+# FITS bitpix. Notice that these are null ops for IEEE floating point format.
+
+procedure wft_set_blank (fits_bitpix, blank, blank_str)
+
+int	fits_bitpix	# the requested FITS bits per pixel
+long	blank		# the FITS integer value of a blank pixel
+char	blank_str[ARB]	# the encoded FITS integer value of a blank pixel
+
+begin
+	switch (fits_bitpix) {
+	case FITS_BYTE:
+	    blank = long (BYTE_BLANK)
+	    call strcpy ("0", blank_str, LEN_BLANK)
+	case FITS_SHORT:
+	    blank = long (SHORT_BLANK)
+	    call strcpy ("-32768", blank_str, LEN_BLANK)
+	case FITS_LONG:
+	    blank = long (LONG_BLANK)
+	    call strcpy ("-2147483648", blank_str, LEN_BLANK)
+	case FITS_REAL:
+	    blank = INDEFL
+	    call strcpy ("", blank_str, LEN_BLANK)
+	case FITS_DOUBLE:
+	    blank = INDEFL
+	    call strcpy ("", blank_str, LEN_BLANK)
+	default:
+	    call error (5, "SET_BLANK: Unknown FITS type.")
+	}
+end
+
+
+# WFT_INIT_CARD_ENCODE -- This procedure initializes the card encoding
+# procedure.  The cards counters are initialized and the number of history cards
+# calculated.
+
+int procedure wft_init_card_encode (im, fits)
+
+# both entry points
+pointer	im	# pointer to the IRAF image
+pointer	fits	# pointer to the WFITS structure
+
+# entry wft_card_encode
+int	wft_card_encode		# entry point
+char	card[LEN_CARD+1]	# string containing the card image
+
+int	cardno, axisno, optiono, hist_ptr, unknown_ptr
+int	nstandard, noptions, stat
+int	wft_standard_card(), wft_option_card(), wft_last_card()
+int	wft_history_card(), wft_unknown_card()
+errchk	wft_standard_card, wft_option_card, wft_history_card
+errchk	wft_unknown_card, wft_last_card
+
+begin
+	# Initialize the card pointers.
+	cardno = 1
+	axisno = 1
+	optiono = 1
+	unknown_ptr = 1
+	hist_ptr = 1
+
+	# Initilaize the card counters.
+	nstandard = 3 + NAXIS(im)
+	noptions = NOPTIONS + nstandard
+
+	return (YES)
+
+
+# WFT_CARD_ENCODE -- Procedure to encode the FITS header parameters into
+# FITS card images.
+
+entry	wft_card_encode (im, fits, card)
+
+	# Fetch the appropriate FITS header card image.
+	if (cardno <= nstandard) {
+	    stat = wft_standard_card (cardno, im, fits, axisno, card)
+	} else if (cardno <= noptions) {
+	    stat = wft_option_card (im, fits, optiono, card)
+	} else if (wft_unknown_card (im, unknown_ptr, card) == YES) {
+	    stat = YES
+	} else if (wft_history_card (im, hist_ptr, card) == YES) {
+	    stat = YES
+	} else {
+	    stat = wft_last_card (card)
+	}
+
+	cardno = cardno + 1
+
+	return (stat)
+end
+
+
+# WFT_TRIMSTR -- Procedure to trim trailing blanks from a fixed size string.
+
+procedure wft_trimstr (instr, outstr, nchars)
+
+char	instr[ARB]	# input string
+char	outstr[ARB]	# output string
+int	nchars		# last character of instr
+
+int	ip
+
+begin
+	call strcpy (instr, outstr, nchars)
+	ip = nchars
+	while (outstr[ip] == ' ')
+	    ip = ip - 1
+	outstr[ip+1] = EOS
+end
diff --git a/pkg/obsolete/fits/fits_wimage.x b/pkg/obsolete/fits/fits_wimage.x
new file mode 100644
index 00000000..539114ea
--- /dev/null
+++ b/pkg/obsolete/fits/fits_wimage.x
@@ -0,0 +1,416 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include <error.h>
+include <imhdr.h>
+include "wfits.h"
+
+# WFT_WRITE_IMAGE -- Procedure to convert IRAF image data to FITS format line by
+# line.
+
+procedure wft_write_image (im, fits, fits_fd)
+
+pointer	im			# IRAF image descriptor
+pointer	fits			# FITS data structure
+int	fits_fd			# FITS file descriptor
+
+int	npix, nlines, npix_record, i, stat, nrecords
+long	v[IM_MAXDIM]
+pointer	tempbuf, buf
+
+int	wft_get_image_line()
+errchk	malloc, mfree, wft_get_image_line, wft_lscale_line, wft_long_line
+errchk	wft_init_write_pixels, wft_write_pixels, wft_write_last_record
+errchk	wft_rscale_line, wft_dscale_line
+
+include "wfits.com"
+
+begin
+	if (NAXIS(im) == 0)
+	    return
+
+	# Initialize.
+	npix = NAXISN(im,1)
+	nlines = 1
+	do i = 2, NAXIS(im)
+	    nlines = nlines * NAXISN(im, i)
+	npix_record = len_record * FITS_BYTE / abs (FITS_BITPIX(fits))
+
+	call amovkl (long(1), v, IM_MAXDIM)
+	switch (FITS_BITPIX(fits)) {
+	case FITS_REAL:
+
+	    # Allocate temporary space.
+	    call malloc (tempbuf, npix, TY_REAL)
+
+	    # Initialize the pixel write.
+	    call wft_init_write_pixels (npix_record, TY_REAL,
+	        FITS_BITPIX(fits), blkfac)
+
+	    # For the time being explicitly turn off ieee NaN mapping.
+	    call ieemapr (NO, NO)
+
+	    # Scale the lines, deal with the blanks via the ieee code which
+	    # is currently turned off, and write the output records.
+
+	    do i = 1, nlines {
+	        iferr (stat =  wft_get_image_line (im, buf, v, PIXTYPE(im))) {
+		    call erract (EA_WARN)
+		    call error (10, "WRT_IMAGE: Error writing IRAF image.")
+	        }
+	        if (stat == EOF )
+		    return
+	        if (stat != npix)
+		    call error (10, "WRT_IMAGE: Error writing IRAF image.")
+	        if (SCALE(fits) == YES)
+	            call wft_rscale_line (buf, Memr[tempbuf], npix,
+		        1. / BSCALE(fits), -BZERO(fits), PIXTYPE(im))
+		else
+		    call wft_real_line (buf, Memr[tempbuf], npix, PIXTYPE(im)) 
+	        call wft_write_pixels (fits_fd, Memr[tempbuf], npix)
+	    }
+
+	    # Free space.
+	    call mfree (tempbuf, TY_REAL)
+
+	case FITS_DOUBLE:
+
+	    # Allocate temporary space.
+	    call malloc (tempbuf, npix, TY_DOUBLE)
+
+	    # Initialize the pixel write.
+	    call wft_init_write_pixels (npix_record, TY_DOUBLE,
+	        FITS_BITPIX(fits), blkfac)
+
+	    # For the time being explicitly turn off ieee NaN mapping.
+	    call ieemapd (NO, NO)
+
+	    # Scale the lines, deal with the blanks via the ieee code which
+	    # is currently turned off, and write the output records.
+
+	    do i = 1, nlines {
+	        iferr (stat =  wft_get_image_line (im, buf, v, PIXTYPE(im))) {
+		    call erract (EA_WARN)
+		    call error (10, "WRT_IMAGE: Error writing IRAF image.")
+	        }
+	        if (stat == EOF )
+		    return
+	        if (stat != npix)
+		    call error (10, "WRT_IMAGE: Error writing IRAF image.")
+	        if (SCALE(fits) == YES)
+	            call wft_dscale_line (buf, Memd[tempbuf], npix,
+		        1. / BSCALE(fits), -BZERO(fits), PIXTYPE(im))
+		else
+		    call wft_double_line (buf, Memd[tempbuf], npix,
+		        PIXTYPE(im)) 
+	        call wft_write_pixels (fits_fd, Memd[tempbuf], npix)
+	    }
+
+	    # Free space.
+	    call mfree (tempbuf, TY_DOUBLE)
+
+	default:
+
+	    # Allocate temporary space.
+	    call malloc (tempbuf, npix, TY_LONG)
+
+	    # Scale the line, deal with the blanks, and write the output
+	    # record. At the moement blanks are not dealt with.
+
+	    call wft_init_write_pixels (npix_record, TY_LONG, FITS_BITPIX(fits),
+	         blkfac)
+	    do i = 1, nlines {
+	        iferr (stat =  wft_get_image_line (im, buf, v, PIXTYPE(im))) {
+		    call erract (EA_WARN)
+		    call error (10, "WRT_IMAGE: Error writing IRAF image.")
+	        }
+	        if (stat == EOF )
+		    return
+	        if (stat != npix)
+		    call error (10, "WRT_IMAGE: Error writing IRAF image.")
+	        if (SCALE(fits) == YES)
+	            call wft_lscale_line (buf, Meml[tempbuf], npix,
+		        1. / BSCALE(fits), -BZERO(fits), PIXTYPE(im))
+	        else
+		    call wft_long_line (buf, Meml[tempbuf], npix, PIXTYPE(im)) 
+	        # call map_blanks (im, Meml[tempbuf], blank)
+	        call wft_write_pixels (fits_fd, Meml[tempbuf], npix)
+	    }
+	    # Free space.
+	    call mfree (tempbuf, TY_LONG)
+	}
+
+	# Write the final record.
+	call wft_write_last_record (fits_fd, nrecords)
+	if (short_header == YES || long_header == YES) {
+	    call printf ("%d  Data logical (2880 byte) records written\n")
+	        call pargi (nrecords)
+	}
+end
+
+
+# WFT_GET_IMAGE_LINE -- Procedure to fetch the next image line.
+
+int procedure wft_get_image_line (im, buf, v, datatype)
+
+pointer	im			# IRAF image descriptor
+pointer	buf			# pointer to image line
+long	v[ARB]			# imio dimension descriptor
+int	datatype		# IRAF image data type
+
+int	npix
+int	imgnll(), imgnlr(), imgnld(), imgnlx()
+errchk	imgnll, imgnlr, imgnld, imgnlx
+
+begin
+	switch (datatype) {
+	case TY_SHORT, TY_INT, TY_USHORT, TY_LONG:
+	    npix = imgnll (im, buf, v)
+	case TY_REAL:
+	    npix = imgnlr (im, buf, v)
+	case TY_DOUBLE:
+	    npix = imgnld (im, buf, v)
+	case TY_COMPLEX:
+	    npix = imgnlx (im, buf, v)
+	default:
+	    call error (11, "GET_IMAGE_LINE: Unknown IRAF image type.")
+	}
+
+	return (npix)
+end
+
+
+# WFT_RSCALE_LINE -- This procedure converts the IRAF data to type real
+# and scales by the FITS parameters bscale and bzero.
+
+procedure wft_rscale_line (buf, outbuffer, npix, bscale, bzero, datatype)
+
+pointer	buf			# pointer to IRAF image line
+real	outbuffer[ARB]		# FITS integer buffer
+int	npix			# number of pixels
+double	bscale, bzero		# FITS bscale and bzero parameters
+int	datatype		# data type of image
+
+errchk	achtlr, altadr, amovr, achtdr, acthxr
+
+begin
+	switch (datatype) {
+	case TY_SHORT, TY_INT, TY_LONG, TY_USHORT:
+	    call achtlr (Meml[buf], outbuffer, npix)
+	    call altadr (outbuffer, outbuffer, npix, bzero, bscale)
+	case TY_REAL:
+	    call amovr (Memr[buf], outbuffer, npix)
+	    call altadr (outbuffer, outbuffer, npix, bzero, bscale)
+	case TY_DOUBLE:
+	    call achtdr (Memd[buf], outbuffer, npix)
+	    call altadr (outbuffer, outbuffer, npix, bzero, bscale)
+	case TY_COMPLEX:
+	    call achtxr (Memx[buf], outbuffer, npix)
+	    call altadr (outbuffer, outbuffer, npix, bzero, bscale)
+	default:
+	    call error (12, "WFT_RSCALE_LINE: Unknown IRAF image type.")
+	}
+end
+
+
+# WFT_DSCALE_LINE -- This procedure converts the IRAF data to type double with
+# after scaling by the FITS parameters bscale and bzero.
+
+procedure wft_dscale_line (buf, outbuffer, npix, bscale, bzero, datatype)
+
+pointer	buf			# pointer to IRAF image line
+double	outbuffer[ARB]		# FITS integer buffer
+int	npix			# number of pixels
+double	bscale, bzero		# FITS bscale and bzero parameters
+int	datatype		# data type of image
+
+errchk	achtld, altad, amovd, achtrd, achtxd
+
+begin
+	switch (datatype) {
+	case TY_SHORT, TY_INT, TY_LONG, TY_USHORT:
+	    call achtld (Meml[buf], outbuffer, npix)
+	    call altad (outbuffer, outbuffer, npix, bzero, bscale)
+	case TY_REAL:
+	    call achtrd (Memr[buf], outbuffer, npix)
+	    call altad (outbuffer, outbuffer, npix, bzero, bscale)
+	case TY_DOUBLE:
+	    call amovd (Memd[buf], outbuffer, npix)
+	    call altad (outbuffer, outbuffer, npix, bzero, bscale)
+	case TY_COMPLEX:
+	    call achtxd (Memx[buf], outbuffer, npix)
+	    call altad (outbuffer, outbuffer, npix, bzero, bscale)
+	default:
+	    call error (12, "WFT_DSCALE_LINE: Unknown IRAF image type.")
+	}
+end
+
+
+# WFT_LSCALE_LINE -- This procedure converts the IRAF data to type long with
+# after scaling by the FITS parameters bscale and bzero.
+
+procedure wft_lscale_line (buf, outbuffer, npix, bscale, bzero, datatype)
+
+pointer	buf			# pointer to IRAF image line
+long	outbuffer[ARB]		# FITS integer buffer
+int	npix			# number of pixels
+double	bscale, bzero		# FITS bscale and bzero parameters
+int	datatype		# data type of image
+
+errchk  altal, amovl, altadr, achtrl, altad, achtdl, altax, achtxl
+
+begin
+	switch (datatype) {
+	case TY_SHORT, TY_INT, TY_LONG, TY_USHORT:
+	    call altal (Meml[buf], Meml[buf], npix, bzero, bscale)
+	    call amovl (Meml[buf], outbuffer, npix)
+	case TY_REAL:
+	    call altarl (Memr[buf], outbuffer, npix, bzero, bscale)
+	    #call altadr (Memr[buf], Memr[buf], npix, bzero, bscale)
+	    #call achtrl (Memr[buf], outbuffer, npix)
+	case TY_DOUBLE:
+	    call altad (Memd[buf], Memd[buf], npix, bzero, bscale)
+	    call achtdl (Memd[buf], outbuffer, npix)
+	case TY_COMPLEX:
+	    call altadx (Memx[buf], Memx[buf], npix, bzero, bscale)
+	    call achtxl (Memx[buf], outbuffer, npix)
+	default:
+	    call error (12, "WFT_LSCALE_LINE: Unknown IRAF image type.")
+	}
+end
+
+
+# WFT_REAL_LINE -- This procedure converts the IRAF image line to type long with
+# no scaling.
+
+procedure wft_real_line (buf, outbuffer, npix, datatype)
+
+pointer	buf			# pointer to IRAF image line
+real	outbuffer[ARB]		# buffer of FITS integers
+int	npix			# number of pixels
+int	datatype		# IRAF image datatype
+
+errchk	achtlr, achtdr, amovr, achtxr
+
+begin
+	switch (datatype) {
+	case TY_SHORT, TY_INT, TY_LONG, TY_USHORT:
+	    call achtlr (Meml[buf], outbuffer, npix)
+	case TY_REAL:
+	    call amovr (Memr[buf], outbuffer, npix)
+	case TY_DOUBLE:
+	    call achtdr (Memd[buf], outbuffer, npix)
+	case TY_COMPLEX:
+	    call achtxr (Memx[buf], outbuffer, npix)
+	default:
+	    call error (13, "WFT_REAL_LINE: Unknown IRAF data type.")
+	}
+end
+
+
+# WFT_DOUBLE_LINE -- This procedure converts the IRAF image line to type long
+# with no scaling.
+
+procedure wft_double_line (buf, outbuffer, npix, datatype)
+
+pointer	buf			# pointer to IRAF image line
+double  outbuffer[ARB]		# buffer of FITS integers
+int	npix			# number of pixels
+int	datatype		# IRAF image datatype
+
+errchk	achtld, achtrd, amovd, achtxd
+
+begin
+	switch (datatype) {
+	case TY_SHORT, TY_INT, TY_LONG, TY_USHORT:
+	    call achtld (Meml[buf], outbuffer, npix)
+	case TY_REAL:
+	    call achtrd (Memr[buf], outbuffer, npix)
+	case TY_DOUBLE:
+	    call amovd (Memd[buf], outbuffer, npix)
+	case TY_COMPLEX:
+	    call achtxd (Memx[buf], outbuffer, npix)
+	default:
+	    call error (13, "WFT_DOUBLE_LINE: Unknown IRAF data type.")
+	}
+end
+
+
+# WFT_LONG_LINE -- This procedure converts the IRAF image line to type long with
+# no scaling.
+
+procedure wft_long_line (buf, outbuffer, npix, datatype)
+
+pointer	buf			# pointer to IRAF image line
+long	outbuffer[ARB]		# buffer of FITS integers
+int	npix			# number of pixels
+int	datatype		# IRAF image datatype
+
+errchk	amovl, achtrl, achtdl, achtxl
+
+begin
+	switch (datatype) {
+	case TY_SHORT, TY_INT, TY_LONG, TY_USHORT:
+	    call amovl (Meml[buf], outbuffer, npix)
+	case TY_REAL:
+	    call achtrl (Memr[buf], outbuffer, npix)
+	case TY_DOUBLE:
+	    call achtdl (Memd[buf], outbuffer, npix)
+	case TY_COMPLEX:
+	    call achtxl (Memx[buf], outbuffer, npix)
+	default:
+	    call error (13, "WFT_LONG_LINE: Unknown IRAF data type.")
+	}
+end
+
+
+# ALTARL -- Procedure to linearly scale a real vector into a long vector
+# using double precision constants to preserve precision.
+
+procedure altarl (a, b, npix, k1, k2)
+
+real	a[ARB]		# input vector
+long	b[ARB]		# output vector
+int	npix		# number of pixels
+double	k1, k2		# scaling factors
+
+int	i
+
+begin
+	do i = 1, npix
+	    b[i] = (a[i] + k1) * k2
+end
+
+
+# ALTADR -- Procedure to linearly scale a real vector in double precision
+
+procedure altadr (a, b, npix, k1, k2)
+
+real	a[ARB]		# input vector
+real	b[ARB]		# output vector
+int	npix		# number of pixels
+double	k1, k2		# scaling factors
+
+int	i
+
+begin
+	do i = 1, npix
+	    b[i] = (a[i] + k1) * k2
+end
+
+
+# ALTADX -- Procedure to linearly scale a complex vector in double precision
+
+procedure altadx (a, b, npix, k1, k2)
+
+complex	a[ARB]		# input vector
+complex	b[ARB]		# output vector
+int	npix		# number of pixels
+double	k1, k2		# scaling factors
+
+int	i
+
+begin
+	do i = 1, npix
+	    b[i] = (a[i] + k1) * k2
+end
+
diff --git a/pkg/obsolete/fits/fits_wpixels.x b/pkg/obsolete/fits/fits_wpixels.x
new file mode 100644
index 00000000..11d4fd0a
--- /dev/null
+++ b/pkg/obsolete/fits/fits_wpixels.x
@@ -0,0 +1,162 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include <mach.h>
+include <fset.h>
+include "wfits.h"
+
+# WFT_INIT_WRITE_PIXELS -- This procedure calculates the input and
+# output buffer sizes based in the spp and mii data types and allocates
+# the required space.
+
+procedure wft_init_write_pixels (npix_record, spp_type, bitpix, blkfac)
+
+int	npix_record		# number of data pixels per record
+int	spp_type		# pixel data type
+int	bitpix			# output bits per pixel
+int	blkfac			# fits blocking factor (0 for disk)
+
+# entry wft_write_pixels, wft_write_last_record
+
+int	fd			# output file descriptor
+char	buffer[1]		# input buffer
+int	npix			# number of pixels in the input buffer
+int	nrecords		# number of FITS records written
+
+char	blank, zero
+int	ty_mii, ty_spp, npix_rec, nch_rec, len_mii, sz_rec, nchars, n, nrec
+int	bf, szblk
+pointer	spp, mii, ip, op
+
+int	sizeof(), miilen(), fstati()
+long	note()
+errchk	malloc, mfree, write, miipak, amovc
+data	mii /NULL/, spp/NULL/
+
+begin
+	# Change input parameters into local variables.
+	ty_mii = bitpix
+	ty_spp = spp_type
+	npix_rec = npix_record
+	nch_rec = npix_rec * sizeof (ty_spp)
+	bf = blkfac
+	blank = ' '
+	zero = 0
+
+	# Compute the size of the mii buffer.
+	len_mii = miilen (npix_rec, ty_mii)
+	sz_rec = len_mii * SZ_INT
+
+	# Allocate space for the buffers.
+	if (spp != NULL)
+	    call mfree (spp, TY_CHAR)
+	call malloc (spp, nch_rec, TY_CHAR)
+	if (mii != NULL)
+	    call mfree (mii, TY_INT)
+	call malloc (mii, len_mii, TY_INT)
+
+	op = 0
+	nrec = 0
+
+	return
+
+# WFT_WRITE_PIXELS -- Wft_wrt_pixels gets an image line and places it in the
+# output buffer. When the output buffer is full the data are packed by the mii
+# routines and written to the specified output.
+
+entry	wft_write_pixels (fd, buffer, npix)
+
+	nchars = npix * sizeof (ty_spp)
+	ip = 0
+
+	repeat {
+
+	    # Fill output buffer.
+	    n = min (nch_rec - op, nchars - ip)
+	    call amovc (buffer[1 + ip], Memc[spp + op], n)
+	    ip = ip + n
+	    op = op + n
+
+	    # Write output record.
+	    if (op == nch_rec) {
+		call miipak (Memc[spp], Memi[mii], npix_rec, ty_spp, ty_mii)
+		iferr (call write (fd, Memi[mii], sz_rec)) {
+	            if (ty_spp == TY_CHAR) {
+		        call printf (" File incomplete: %d logical header")
+			    call pargi (nrec)
+			call printf (" (2880 byte) records written\n")
+	                call error (18,
+			    "WRT_RECORD: Error writing header record.")
+	            } else {
+		        call printf (" File incomplete: %d logical data")
+			    call pargi (nrec)
+			call printf (" (2880 byte) records written\n")
+		        call error (19,
+			    "WRT_RECORD: Error writing data record.")
+		    }
+		}
+
+		nrec = nrec + 1
+		op = 0
+	    }
+
+	} until (ip == nchars)
+
+	return
+
+
+# WFT_WRITE_LAST_RECORD -- Procedure to write the last partially filled record
+# to tape. Fill with blanks if header record otherwise fill with zeros.
+
+entry	wft_write_last_record (fd, nrecords)
+
+	if (op != 0) {
+
+	    # Blank or zero fill the last record.
+	    n = nch_rec - op
+	    if (ty_spp == TY_CHAR)
+		call amovkc (blank, Memc[spp + op], n)
+	    else
+		call amovkc (zero, Memc[spp + op], n)
+
+	    # Write last record.
+	    call miipak (Memc[spp], Memi[mii], npix_rec, ty_spp, ty_mii)
+	    iferr (call write (fd, Memi[mii], sz_rec)) {
+	        if (ty_spp == TY_CHAR) {
+		    call printf ("File incomplete: %d logical header")
+			call pargi (nrec)
+		    call printf (" (2880 byte) records written\n")
+	            call error (18,
+			"WRT_LAST_RECORD: Error writing last header record.")
+	        } else {
+		    call printf ("File incomplete: %d logical data")
+			call pargi (nrec)
+		    call printf (" (2880 byte) records written\n")
+		    call error (19,
+			"WRT_LAST_RECORD: Error writing last data record.")
+		}
+	    }
+
+
+	    nrec = nrec + 1
+
+	    # Pad out the record if the blocking is non-standard.
+	    szblk = fstati (fd, F_BUFSIZE) * SZB_CHAR
+	    if ((bf > 0) && mod (szblk, FITS_RECORD) != 0 &&
+	        (ty_spp != TY_CHAR)) {
+		szblk = szblk / SZB_CHAR
+		n = note (fd) - 1
+		if (mod (n, szblk) == 0)
+		    n = 0
+		else
+		    n = szblk - mod (n, szblk)
+		for (op = 1; op <= n; op = op + nch_rec) {
+		    szblk = min (nch_rec, n - op + 1)
+		    call amovkc (zero, Memc[spp], szblk)
+		    #call write (fd, Memc[spp], szblk)
+		}
+	    }
+
+	}
+
+	nrecords = nrec
+end
diff --git a/pkg/obsolete/fits/fits_write.x b/pkg/obsolete/fits/fits_write.x
new file mode 100644
index 00000000..6887ead0
--- /dev/null
+++ b/pkg/obsolete/fits/fits_write.x
@@ -0,0 +1,156 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include <fset.h>
+include	<error.h>
+include <mach.h>
+include <imhdr.h>
+include "wfits.h"
+
+# WFT_WRITE_FITZ -- Procedure to convert a single IRAF file to a FITS file.
+# If the make_image switch is set the header and pixel files are output
+# to the output destination. If the make_image switch is off the header
+# is printed to the standard output.
+
+procedure wft_write_fitz (iraf_file, fits_file)
+
+char	iraf_file[ARB]		# IRAF file name
+char	fits_file[ARB]		# FITS file name
+
+int	fits_fd, chars_rec, nchars, ip, min_lenuserarea
+pointer	im, sp, fits, envstr
+
+int	mtfile(), mtopen(), open(), fnldir(), envfind(), ctoi()
+pointer	immap()
+errchk	immap, imunmap, open, mtopen, close, smark, salloc, sfree
+errchk	delete, wft_write_header, wft_write_image, wft_data_limits
+
+include "wfits.com"
+
+begin
+	# Allocate memory for program data structure.
+	call smark (sp)
+	call salloc (fits, LEN_FITS, TY_STRUCT)
+	call salloc (envstr, SZ_FNAME, TY_CHAR)
+
+	# Construct the old iraf name by removing the directory
+	# specification.
+
+	call imgcluster (iraf_file, IRAFNAME(fits), SZ_FNAME)
+	nchars = fnldir (IRAFNAME(fits), IRAFNAME(fits), SZ_FNAME)
+	call strcpy (iraf_file[nchars+1], IRAFNAME(fits), SZ_FNAME)
+
+	# Open the input image.
+	if (envfind ("min_lenuserarea", Memc[envstr], SZ_FNAME) > 0) {
+	    ip = 1
+	    if (ctoi (Memc[envstr], ip, min_lenuserarea) <= 0)
+		min_lenuserarea = LEN_USERAREA
+	    else
+		min_lenuserarea = max (LEN_USERAREA, min_lenuserarea)
+	} else
+	    min_lenuserarea = LEN_USERAREA
+	im = immap (iraf_file, READ_ONLY, min_lenuserarea)
+
+	# Open the output file. Check whether the output file is a magtape
+	# device or a binary file. If the output file is magtape check
+	# for a legal blocking factor.
+
+	if (make_image == NO)
+	    call strcpy ("dev$null", fits_file, SZ_FNAME)
+
+	if (mtfile (fits_file) == YES) {
+	    chars_rec = (blkfac * len_record * FITS_BYTE) / (SZB_CHAR *
+	        NBITS_BYTE)
+	    fits_fd = mtopen (fits_file, WRITE_ONLY, chars_rec)
+	} else
+	    fits_fd = open (fits_file, NEW_FILE, BINARY_FILE)
+
+	# Write header and image.
+	iferr {
+
+	    if (short_header == YES || long_header == YES) {
+	        if (make_image == YES) {
+		    call printf (" -> %s ")
+		        call pargstr (fits_file)
+	        }
+		if (long_header == YES)
+		    call printf ("\n")
+	    }
+	    call flush (STDOUT)
+
+	    call wft_write_header (im, fits, fits_fd)
+	    if (make_image == YES)
+		call wft_write_image (im, fits, fits_fd)
+
+	    if (long_header == YES)
+	        call printf ("\n")
+
+	} then {
+
+	    # Print the error message.
+	    call flush (STDOUT)
+	    call erract (EA_WARN)
+
+	    # Close files and cleanup.
+	    call imunmap (im)
+	    call close (fits_fd)
+	    if (make_image == NO)
+	        call delete (fits_file)
+	    call sfree (sp)
+
+	    # Assert an error.
+	    call erract (EA_ERROR)
+
+	} else {
+
+	    # Close files and cleanup.
+	    call imunmap (im)
+	    call close (fits_fd)
+	    if (make_image == NO)
+	        call delete (fits_file)
+	    call sfree (sp)
+	}
+
+end
+
+
+# WFT_DATA_LIMITS -- Procedure to calculate the maximum and minimum data values
+# in an IRAF image. Values are only calculated if the max and min are unknown
+# or the image has been modified since the last values were calculated.
+
+procedure wft_data_limits (im, irafmin, irafmax)
+
+pointer	im		# image pointer
+real	irafmin		# minimum picture value
+real	irafmax		# maximum picture value
+
+int	npix
+long	v[IM_MAXDIM]
+pointer	buf
+real	maxval, minval
+int	imgnlr()
+errchk	imgnlr
+
+begin
+	# Compute the data minimum and maximum if the image values
+	# are undefined out-of-date.
+
+	if (LIMTIME(im) < MTIME(im) && NAXIS(im) > 0) {
+
+	    irafmax = -MAX_REAL
+	    irafmin = MAX_REAL
+	    npix = NAXISN(im,1)
+
+	    call amovkl (long(1), v, IM_MAXDIM)
+	    while (imgnlr (im, buf, v) != EOF) {
+	        call alimr (Memr[buf], npix, minval, maxval)
+	        irafmin = min (irafmin, minval)
+	        irafmax = max (irafmax, maxval)
+	    }
+
+	} else {
+
+	    irafmax = IM_MAX(im)
+	    irafmin = IM_MIN(im)
+
+	}
+end
diff --git a/pkg/obsolete/fits/mkpkg b/pkg/obsolete/fits/mkpkg
new file mode 100644
index 00000000..dd583c01
--- /dev/null
+++ b/pkg/obsolete/fits/mkpkg
@@ -0,0 +1,23 @@
+# Fits Library
+
+$checkout libpkg.a ../
+$update   libpkg.a
+$checkin  libpkg.a ../
+$exit
+
+libpkg.a:
+	fits_cards.x	wfits.com wfits.h <imhdr.h>
+	fits_params.x	wfits.h <time.h>
+	fits_read.x	rfits.com rfits.h <error.h> <fset.h> <imhdr.h>
+	fits_rheader.x	rfits.com rfits.h rfits.com <ctype.h> <imhdr.h>\
+			<imio.h> <mach.h>
+	fits_rimage.x	rfits.com rfits.h <imhdr.h> <fset.h> <mach.h>
+	fits_rpixels.x	<fset.h> <mach.h> <mii.h>
+	fits_wheader.x	wfits.com wfits.h <fset.h> <imhdr.h> <mach.h>
+	fits_wimage.x	wfits.com wfits.h <error.h> <imhdr.h>
+	fits_wpixels.x	wfits.h <fset.h> <mach.h>
+	fits_write.x	<error.h> wfits.com wfits.h <fset.h> <imhdr.h> <mach.h>
+	ranges.x	<mach.h> <ctype.h>
+	t_rfits.x	rfits.com rfits.h <error.h> <fset.h>
+	t_wfits.x	wfits.com wfits.h <error.h> <fset.h> <mach.h>
+	;
diff --git a/pkg/obsolete/fits/ranges.x b/pkg/obsolete/fits/ranges.x
new file mode 100644
index 00000000..b3812cd1
--- /dev/null
+++ b/pkg/obsolete/fits/ranges.x
@@ -0,0 +1,234 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<mach.h>
+include	<ctype.h>
+
+define	FIRST	1		# Default starting range
+define	LAST	MAX_INT		# Default ending range
+define	STEP	1		# Default step
+
+# DECODE_RANGES -- Parse a string containing a list of integer numbers or
+# ranges, delimited by either spaces or commas.  Return as output a list
+# of ranges defining a list of numbers, and the count of list numbers.
+# Range limits must be positive nonnegative integers.  ERR is returned as
+# the function value if a conversion error occurs.  The list of ranges is
+# delimited by a single NULL.
+
+int procedure decode_ranges (range_string, ranges, max_ranges, nvalues)
+
+char	range_string[ARB]	# Range string to be decoded
+int	ranges[3, max_ranges]	# Range array
+int	max_ranges		# Maximum number of ranges
+int	nvalues			# The number of values in the ranges
+
+int	ip, nrange, first, last, step, ctoi()
+
+begin
+	ip = 1
+	nvalues = 0
+
+	do nrange = 1, max_ranges - 1 {
+	    # Defaults to all positive integers
+	    first = FIRST
+	    last = LAST
+	    step = STEP
+
+	    # Skip delimiters
+	    while (IS_WHITE(range_string[ip]) || range_string[ip] == ',')
+		ip = ip + 1
+
+	    # Get first limit.
+	    # Must be a number, '-', 'x', or EOS.  If not return ERR.
+	    if (range_string[ip] == EOS) {			# end of list
+		if (nrange == 1) {
+		    # Null string defaults
+		    ranges[1, 1] = first
+		    ranges[2, 1] = last
+		    ranges[3, 1] = step
+		    ranges[1, 2] = NULL
+	    	    nvalues = nvalues + abs (last-first) / step + 1
+		    return (OK)
+		} else {
+		    ranges[1, nrange] = NULL
+		    return (OK)
+		}
+	    } else if (range_string[ip] == '-')
+		;
+	    else if (range_string[ip] == 'x')
+		;
+	    else if (IS_DIGIT(range_string[ip])) {		# ,n..
+		if (ctoi (range_string, ip, first) == 0)
+		    return (ERR)
+	    } else
+		return (ERR)
+
+	    # Skip delimiters
+	    while (IS_WHITE(range_string[ip]) || range_string[ip] == ',')
+		ip = ip + 1
+
+	    # Get last limit
+	    # Must be '-', or 'x' otherwise last = first.
+	    if (range_string[ip] == 'x')
+		;
+	    else if (range_string[ip] == '-') {
+		ip = ip + 1
+	        while (IS_WHITE(range_string[ip]) || range_string[ip] == ',')
+		    ip = ip + 1
+		if (range_string[ip] == EOS)
+		    ;
+		else if (IS_DIGIT(range_string[ip])) {
+		    if (ctoi (range_string, ip, last) == 0)
+		        return (ERR)
+		} else if (range_string[ip] == 'x')
+		    ;
+		else
+		    return (ERR)
+	    } else
+		last = first
+
+	    # Skip delimiters
+	    while (IS_WHITE(range_string[ip]) || range_string[ip] == ',')
+		ip = ip + 1
+
+	    # Get step.
+	    # Must be 'x' or assume default step.
+	    if (range_string[ip] == 'x') {
+		ip = ip + 1
+	        while (IS_WHITE(range_string[ip]) || range_string[ip] == ',')
+		    ip = ip + 1
+		if (range_string[ip] == EOS)
+		    ;
+		else if (IS_DIGIT(range_string[ip])) {
+		    if (ctoi (range_string, ip, step) == 0)
+		        ;
+		} else if (range_string[ip] == '-')
+		    ;
+		else
+		    return (ERR)
+	    }
+
+	    # Output the range triple.
+	    ranges[1, nrange] = first
+	    ranges[2, nrange] = last
+	    ranges[3, nrange] = step
+	    nvalues = nvalues + abs (last-first) / step + 1
+	}
+
+	return (ERR)					# ran out of space
+end
+
+
+# GET_NEXT_NUMBER -- Given a list of ranges and the current file number,
+# find and return the next file number.  Selection is done in such a way
+# that list numbers are always returned in monotonically increasing order,
+# regardless of the order in which the ranges are given.  Duplicate entries
+# are ignored.  EOF is returned at the end of the list.
+
+int procedure get_next_number (ranges, number)
+
+int	ranges[ARB]		# Range array
+int	number			# Both input and output parameter
+
+int	ip, first, last, step, next_number, remainder
+
+begin
+	# If number+1 is anywhere in the list, that is the next number,
+	# otherwise the next number is the smallest number in the list which
+	# is greater than number+1.
+
+	number = number + 1
+	next_number = MAX_INT
+
+	for (ip=1;  ranges[ip] != NULL;  ip=ip+3) {
+	    first = min (ranges[ip], ranges[ip+1])
+	    last = max (ranges[ip], ranges[ip+1])
+	    step = ranges[ip+2]
+	    if (number >= first && number <= last) {
+		remainder = mod (number - first, step)
+		if (remainder == 0)
+		    return (number)
+		if (number - remainder + step <= last)
+		    next_number = number - remainder + step
+	    } else if (first > number)
+		next_number = min (next_number, first)
+	}
+
+	if (next_number == MAX_INT)
+	    return (EOF)
+	else {
+	    number = next_number
+	    return (number)
+	}
+end
+
+
+# GET_PREVIOUS_NUMBER -- Given a list of ranges and the current file number,
+# find and return the previous file number.  Selection is done in such a way
+# that list numbers are always returned in monotonically decreasing order,
+# regardless of the order in which the ranges are given.  Duplicate entries
+# are ignored.  EOF is returned at the end of the list.
+
+int procedure get_previous_number (ranges, number)
+
+int	ranges[ARB]		# Range array
+int	number			# Both input and output parameter
+
+int	ip, first, last, step, next_number, remainder
+
+begin
+	# If number-1 is anywhere in the list, that is the previous number,
+	# otherwise the previous number is the largest number in the list which
+	# is less than number-1.
+
+	number = number - 1
+	next_number = 0
+
+	for (ip=1;  ranges[ip] != NULL;  ip=ip+3) {
+	    first = min (ranges[ip], ranges[ip+1])
+	    last = max (ranges[ip], ranges[ip+1])
+	    step = ranges[ip+2]
+	    if (number >= first && number <= last) {
+		remainder = mod (number - first, step)
+		if (remainder == 0)
+		    return (number)
+		if (number - remainder >= first)
+		    next_number = number - remainder
+	    } else if (last < number) {
+		remainder = mod (last - first, step)
+		if (remainder == 0)
+		    next_number = max (next_number, last)
+		else if (last - remainder >= first)
+		    next_number = max (next_number, last - remainder)
+	    }
+	}
+
+	if (next_number == 0)
+	    return (EOF)
+	else {
+	    number = next_number
+	    return (number)
+	}
+end
+
+
+# IS_IN_RANGE -- Test number to see if it is in range.
+
+bool procedure is_in_range (ranges, number)
+
+int	ranges[ARB]		# Range array
+int	number			# Number to be tested against ranges
+
+int	ip, first, last, step
+
+begin
+	for (ip=1;  ranges[ip] != NULL;  ip=ip+3) {
+	    first = min (ranges[ip], ranges[ip+1])
+	    last = max (ranges[ip], ranges[ip+1])
+	    step = ranges[ip+2]
+	    if (number >= first && number <= last)
+		if (mod (number - first, step) == 0)
+		    return (true)
+	}
+
+	return (false)
+end
diff --git a/pkg/obsolete/fits/rfits.com b/pkg/obsolete/fits/rfits.com
new file mode 100644
index 00000000..08f44c0e
--- /dev/null
+++ b/pkg/obsolete/fits/rfits.com
@@ -0,0 +1,18 @@
+
+# FITS reader common
+
+int	len_record		# Length of FITS records in bytes
+int	data_type		# Output data type
+real	blank			# Blank value
+real	fe			# Maximum size in megabytes for scan mode
+
+# Option flags
+
+int	make_image		# Create an IRAF image
+int	long_header		# Print a long header (FITS header cards)
+int	short_header		# Print a short header (Title and size)
+int	scale			# Scale the data
+int	old_name		# Use old IRAF name?
+
+common	/rfitscom/ len_record, data_type, blank, fe, make_image, long_header,
+	short_header, scale, old_name
diff --git a/pkg/obsolete/fits/rfits.h b/pkg/obsolete/fits/rfits.h
new file mode 100644
index 00000000..8afc5b0b
--- /dev/null
+++ b/pkg/obsolete/fits/rfits.h
@@ -0,0 +1,80 @@
+# FITS Definitions
+
+# The FITS standard readable by the FITS reader using these definitions:
+# 
+# 1.  8 bits / byte
+# 2.  ASCII character code
+# 3.  MII integer data format (i.e. 8 bit unsigned integers and 16 and 32
+#     bit signed twos complement integers with most significant bytes first.)
+# 4.  IEEE 32 and 64 bit floating point format
+#
+#
+# The following deviations from the FITS standard are allowed:
+# 
+# 1. The number of FITS bytes per record is normally 2880 or up to 10 times
+#    2880 bytes but may be arbitrarily specified  by the user.
+
+# Define the bits per pixel, precision and byte order of the basic FITS types
+
+define	FITS_RECORD	2880	# number of bytes in a standard FITS record
+
+define	FITS_BYTE	8	# Bits in a FITS byte
+define	FITS_SHORT	16	# Bits in a FITS short
+define	FITS_LONG	32	# Bits in a FITS long
+define	FITS_REAL	-32	# Bits in a FITS real * -1
+define	FITS_DOUBLE	-64	# Bits in a FITS double * -1
+
+define	FITSB_PREC	3	# Decimal digits of precision in a FITS byte
+define	FITSS_PREC	5	# Decimal digits of precision in a FITS short
+define	FITSL_PREC	10	# Decimal digits of precision in a FITS long
+
+define	LSBF		NO	# Least significant byte first
+
+# Define the basic format of a FITS cardimage
+
+define	LEN_CARD	80	# Length of FITS card in characters
+define	COL_VALUE	11	# Starting column for parameter values
+
+
+# FITS standards not recognized currently by IRAF.
+#
+# 1. SIMPLE	SIMPLE = 'F' not implemented, file skipped
+# 2. GROUPS	Group data not currently implemented, file skippped
+
+
+# Values for the following quantities are stored in the structure below.
+
+define	LEN_FITS	(15 + SZ_FNAME + 1)
+
+define	FITS_BSCALE	Memd[P2D($1)]	  # FITS scaling parameter
+define	FITS_BZERO	Memd[P2D($1+2)]   # FITS zero parameter
+define	BLANK_VALUE	Meml[P2L($1+4)]   # Blank value
+define	BLANKS		Memi[$1+5]	  # YES if blank keyword in header
+define	BITPIX		Memi[$1+6]	  # Bits per pixel (Must be an MII type)
+define	SCALE		Memi[$1+7]	  # Scale the data?
+define	SIMPLE		Memi[$1+8]	  # Standard FITS format 
+define	NRECORDS	Memi[$1+9]	  # Number of FITS logical records
+define	IRAFNAME	Memc[P2C($1+12)]  # Old IRAF name
+
+# Mapping of additional IRAF header parameters
+
+define	PIXTYPE		IM_PIXTYPE($1)
+define	NBPIX		IM_NBPIX($1)
+define	IRAFMAX		IM_MAX($1)
+define	IRAFMIN		IM_MIN($1)
+define	LIMTIME		IM_LIMTIME($1)
+define	LEN_USERAREA	28800
+
+# Mapping of FITS Keywords to IRAF image header
+
+define	NAXIS	 IM_NDIM($1)
+define	NAXISN	 IM_LEN($1,$2)
+define	OBJECT	 IM_TITLE($1)
+define	HISTORY	 IM_HISTORY($1)
+define	UNKNOWN	 Memc[($1+IMU-1)*SZ_STRUCT+1]	# All unrecognized keywords
+						# are stored here
+# Miscellaneous definitions.
+
+define	SZ_OBJECT	SZ_IMTITLE
+define	SZ_HISTORY	SZ_IMHIST
+define	SZ_FCTYPE	SZ_CTYPE
diff --git a/pkg/obsolete/fits/structure.hlp b/pkg/obsolete/fits/structure.hlp
new file mode 100644
index 00000000..715ef185
--- /dev/null
+++ b/pkg/obsolete/fits/structure.hlp
@@ -0,0 +1,363 @@
+.help fits Apr84 "Program Structure"
+.sh
+RFITS Structure Chart
+
+.nf
+t_rfits ()
+# Returns when file list is satisfied or if EOT is encountered
+# Errors from called routines are trapped and printed as a warning.
+
+	read_fits (fitsfile, iraffile)
+	# Returns OK or EOF
+
+		read_header (fits_fd, fits, im)
+		# Returns OK or EOF
+
+	    		decode_fits_card (fits, im, card)
+			# Returns YES or NO
+
+				get_fits_string (card, str, maxchar)
+
+		read_image (fits_fd, fits, im)
+		# Invokes error handler if EOF is encountered
+
+			set_image_header (fits, im)
+
+			set_coords (im)
+			# sets the coordinate transformations to a 1 to 1
+			# mapping
+
+			init_read_pixels (npix_record, bitpix, lsbf, spp_type)
+			# Returns OK
+
+			put_imageline (im, bufptr, v, pixel_type)
+			# Outputs line to image
+
+			read_pixels (fd, buffer, npix)
+			# Returns EOF or the number of pixels converted
+
+			map_blanks (a, blanks, im)
+
+			scale_line (line, bufptr, npix, bscale, bzero,
+				    pixel_type)
+			# Converts the pixels to the output data type after
+			# applying bscale and bzero to the data
+
+			change_line(line, bufptr, npix, pixel_type)
+			# Changes the FITS integers to the output pixel_type
+			# without scaling.
+.fi
+
+.sh
+RFITS Structure Summary
+
+.ls 4 t_rfits
+The main procedure reads the control parameters.
+The files to be read and converted are calculated from the specified source
+and file list.  A loop through
+the files determines the specific input source names and output filenames
+and calls READ_FITS for each conversion.
+.ls read_fits
+The input source is opened and the output image header file is created.
+If only the FITS header is to be listed then a temporary image header file
+is created.  The FITS header is read and decoded into the IRAF image
+header by READ_HEADER.  If the image is to be read then MAKE_IMAGE is
+called.  Finally, all files are closed.  If a temporary image header file
+was created it is deleted.
+.ls read_header
+Some initialization is done on the IRAF header.
+The FITS header cards are read one at a time.  If EOF is encountered
+then EOF is returned.  If a long header listing has been specified then
+the card is printed.  The card is passed to DECODE_FITS_CARD.  If
+DECODE_FITS_CARD returns YES for the END card then the loop exits.  If a
+short header listing has bee specified then the title and image size
+is printed.  The routine returns OK if the header was successfully
+interpreted or EOF if encountered.  All other errors are returned
+via the error handler.
+.ls decode_fits_card
+A series of STRMATCHes are made against the recognized FITS keywords.
+If a match is found the possible actions are:
+.ls
+Convert a value to the IRAF image header.  The conversions are defined in
+fits.h
+.le
+.ls
+Invoke an error condition
+.le
+.ls
+Return a YES status if the keyword is END
+.le
+
+Otherwise, the card is concatenated to the User Area of the IRAF image
+header.  If the keyword was not END then NO is returned.
+.ls get_fits_string
+The string field for a keyword with a string value is extracted.  The
+leading and trailing quotes are removed as well as trailing blanks.
+The EOS is marked by either ', /, or the end of the card.
+.le
+.le
+.le
+.ls read_image
+The FITS image pixels are converted to an IRAF image file.
+The image file header is set.
+The lines of the image are converted one at a time.
+Each line is checked for blank values.
+When all the image data has been read the number of blank pixels encounter
+is printed (unless the value is zero).
+.ls set_image_header
+The pixel type for the IRAF image is set to the user specified type.
+If no type has been specified then the type is determined from the number
+of bits per pixel given in the FITS header.
+.le
+.ls set_coords
+Sets the coordinate transformation parameters to a 1 to 1 transformation.
+.le
+.ls init_read_pixels
+The pixel reading routine is initialized.  The parameters are the number
+of pixels per record, the number of bits per pixel which must be a
+valid MII type, a byte order flag, and the SPP data type to be converted
+to.  In the FITS reader the byte order is specified to be most significant
+byte first and the SPP data type is TY_LONG.
+.le
+.ls put_imageline
+Put_imageline outputs a single line of the FITS image to the IRAF image.
+.le
+.ls read_pixels
+The pixels are read into a record buffer.  Data conversion is accomplished
+with the MII routines since FITS format is the same as MII format.  The
+specified number of pixels is returned in the specified buffer.
+.le
+.ls map_blanks
+Pixels having the blank value as given in the FITS header are added into
+the bad pixel count in the image header.  This routine will someday handle
+mapping of bad pixels in a more detailed manner.
+.le
+.ls scale_line
+The FITS integers from tape are converted to the output IRAF data type by
+applying the FITS scaling parameters BSCALE and BZERO.
+.le
+.ls change_line
+The FITS integers are converted directly to the output IRAF data type
+without scaling (BSCALE = 1.0 and BZERO = 0.0).
+.le
+.le
+.le
+
+.sh
+WFITS Structure Chart
+
+.nf
+t_wfits()
+# Returns when the input file list is satisfied.
+# Errors from called routines are trapped, an error message is issued and
+# wfits terminates.
+
+	wrt_fits (iraf_file, fits_file)
+
+		data_limits (im)
+
+		wrt_header (im, fits, fits_fd)
+
+			set_bitbix (bitpix, pixtype, data_bitpix)
+
+			iraf_type (pixtype, pixstring)
+
+			set_scale (fits_bitpix, datamax, datamin, bscale, bzero)
+
+				tape_limits (fits_bitpix, tapemax, tapemin)
+
+			set_blank (fits_bitpix, blank, blankstr)
+
+			init_wrt_pixels (npix_record, spp_type, mii_type)
+
+			init_card_encode (im, fits)
+			# Returns YES
+
+			card_encode (im, fits, card)
+			# Returns YES if valid FITS card
+
+				get_standard_card (cardno, im, fits, axisno,
+						   card
+				# Returns YES or NO
+
+				get_option_card (im, fits, optiono, card)
+				# Returns YES or NO
+
+				get_coord_card (im, coordno, coordaxis, card)
+				# Returns YES or NO
+
+				get_history_card (im, histptr, card)
+				# Returns YES or NO
+
+				get_end_card (card)
+				# Returns YES or NO
+
+			wrt_pixels (fits_fd, card, len_card)
+
+			wrt_last_record (fits_fd)
+
+		wrt_image (im, fits, fits_fd)
+
+			init_wrt_pixels (npix_record, spp_type, mii_type)
+
+			get_image_line (im, bufptr, v, pixtype)
+			# Returns EOF or number of pixels in a line
+
+			scale_line (bufptr, long_array, npix, bscale, bzero,
+				    pixtype)
+
+			long_line (bufptr, long_array, npix, pixtype)
+
+			map_blanks (im, long_array, blank)
+
+			wrt_pixels (fits_fd, long_array, npix)
+
+			wrt_last_record (fits_fd)
+.fi
+.sh
+WFITS Structure Summary
+
+.ls t_wfits
+The main procedure reads the control parameters. The files to be read and
+converted are calculated from the specified source and file list. A loop
+through the files determines the specific input source names and output
+file names and calls WRT_FITS for each file conversion. Write errors are trapped
+and cause termination of the program.
+.ls wrt_fits
+The input source is opened. If the make_image switch is set the output
+destination is opened. The IRAF image header parameters are encoded into
+the FITS header and both header and pixels are written to the output
+destination.
+DATA_LIMITS is called if the autoscale switch is enabled and the IRAF image
+data maximum and minimum values are unknown or the image has been modified
+since they were last calculated.
+If the make_image switch is turned off the FITS header is printed
+on the standard output and a temporary output file is created. Finally
+all the files are closed. If a temporary file was created it is deleted.
+.ls data_limits
+DATA_LIMITS calculates the minimum and maximum data values in an IRAF image.
+The calculation is made only if these values are
+undefined or the image has been modified since the last values were
+calculated.
+.le
+.ls wrt_header
+Some initialization is done on the FITS header. The appropriate FITS bitpix,
+bzero and bscale factors, and the tape value for blank pixels are calculated
+for each image.
+The FITS header cards are encoded one at a time. If a long_header listing
+has been specified then the FITS card is printed on the standard output.
+If a short_header listing
+is specified then the title and image size only are printed. Encoding terminates
+when the FITS END keyword is encountered. Partial header records are filled
+with blanks before being written to the output destination.
+.ls set_bitpix
+The value of the FITS bitpix is calculated. If the user given bitpix is
+not a valid FITS bitpix, SET_BITPIX uses the precision of the IRAF image data
+to set bitpix.
+.le
+.ls iraf_type
+The IRAF datatype value is set to either INTEGER, FLOATING or COMPLEX.
+.le
+.ls set_scale
+The bscale and bzero values are calculated from the IRAF minimum and maximum
+data values and the FITS bitpix.
+.ls tape_limits
+The maximum and minimum FITS integer values are calculated from the FITS bitpix.
+.le
+.le
+.ls set_blank
+The FITS integer value for a blank pixel is calculated from the FITS bitpix.
+.le
+.ls init_wrt_pixels
+The pixel writing routine is initialized. The parameters are the number of
+pixels per record, the spp_type and the number of bits per pixel which
+must be a valid mii type. For ASCII header card images the number of pixels
+per record is normally 2880 and the bits per pixel is 8. The spp type is
+always TY_CHAR.
+.le
+.ls init_card_encode
+The card encoding procedure is initialized. The number of standard keyword,
+optional keywords, coordinate transformation keywords, and history keywords
+are calculated.
+.le
+.ls card_encode
+An eighty character FITS format string is created for each permitted
+FITS keyword.
+.ls get_standard_card
+The minimum required FITS header parameters, SIMPLE, BITPIX, NAXIS and NAXIS#
+are encode into FITS cards.
+.le
+.ls get_option_card
+A set of optional FITS  parameters are encoded into FITS cards.
+At present the permitted keywords are BSCALE, BZERO, BUNIT, BLANK
+OBJECT, ORIGIN, DATE, IRAFMAX, IRAFMIN, IRAF-B/P, and IRAFTYPE.
+The BLANK card is only encoded if the number of bad pixels in the 
+IRAF image header is nonzero. BUNIT and OBJECT cards are only encoded if
+if the appropriate strings in the IRAF image header are defined.
+.le
+.ls get_coord_card
+The coordinate transformation parameters are encoded into FITS header cards.
+.le
+.ls get_history_card
+The IRAF history string is encoded into FITS header card(s). A maximum of
+seven HISTORY cards are currently permitted by imio.
+.le
+.ls get_end_card
+The FITS end card in encoded.
+.le
+.le
+.ls wrt_pixels
+The FITS card images are read into a record buffer. When the buffer is
+full data conversion is accomplished using the mii routines since FITS
+format is the same as mii format. After data conversion the record
+buffer is written to the output destination.
+.le
+.ls wrt_last_record
+The last partially filled header record is padded with blanks and written to
+the output destination.
+.le
+.le
+.ls wrt_image
+The IRAF image pixels are converted to FITS image format. The lines of the
+IRAF image are converted one at a time. Each line is scaled (if scaling is
+enabled and appropriate) and the IRAF pixels are converted to long integers.
+WRT_PIXELS is called to convert the integers to the appropriate FITS
+output type using the MII routines.
+.ls init_wrt_pixels
+The pixel writing routine is initialized. The parameters are the number of
+pixels per output record, the spp type and the number of FITS bits per pixel
+which must be a valid mii type. The number of pixels per output record is
+2880, 1440 or 720 for FITS bitpix of 8, 16 or 32 respectively. The spp type
+is always TY_LONG.
+.le
+.ls get_image_line
+A single IRAF image line is read into an internal buffer.
+.le
+.ls scale_line
+The IRAF image line is scaled by the FITS BSCALE and
+BZERO scaling factors and converted to long integers.
+SCALE_LINE is called if the scale switch is set.
+.le
+.ls long_line
+The IRAF image data values are converted directly to long integers.
+LONG_LINE is called if scaling switch is turned off.
+.le
+.ls map_blanks
+This function will eventually map IRAF blank pixels to the appropriate
+FITS integer for blanks. Implementation of this function is will occur
+after the imio mods.
+.le
+.ls wrt_pixels
+The scaled IRAF image lines are read into a record buffer. When the buffer is
+full data conversion is accomplished using the MII routines since FITS
+format is the same as mii format. After data conversion the record
+buffer is written to the output destination.
+.le
+.ls wrt_last_record
+The last partially full data record is padded with zeros and written to the
+output destination.
+.le
+.le
+.le
+.le
+.endhelp
diff --git a/pkg/obsolete/fits/t_rfits.x b/pkg/obsolete/fits/t_rfits.x
new file mode 100644
index 00000000..be4449ac
--- /dev/null
+++ b/pkg/obsolete/fits/t_rfits.x
@@ -0,0 +1,184 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include <error.h>
+include <fset.h>
+include "rfits.h"
+
+define	MAX_RANGES	100		# the maxium number of ranges
+define	NTYPES		7		# the number of image data types
+
+# RFITS -- Read FITS format data.  Further documentation given in rfits.hlp
+
+procedure t_rfits()
+
+char	infile[SZ_FNAME]		# fits file
+char	file_list[SZ_LINE]		# list of tape files
+char	outfile[SZ_FNAME]		# IRAF file
+char	in_fname[SZ_FNAME]		# input file name
+char	out_fname[SZ_FNAME]		# output file name
+
+int	len_inlist, len_outlist
+int	range[MAX_RANGES*3+1], file_number, offset, stat
+pointer	inlist, outlist
+
+bool	clgetb()
+char	clgetc()
+int	rft_get_image_type(), clgeti(), mtfile(), strlen(), btoi(), fntlenb()
+int	rft_read_fitz(), decode_ranges(), get_next_number(), fntgfnb(), fstati()
+int	mtneedfileno(), fntrfnb()
+pointer	fntopnb()
+real	clgetr(), rft_fe()
+
+include	"rfits.com"
+
+begin
+	# Set up the standard output to flush on a newline.
+	if (fstati (STDOUT, F_REDIR) == NO)
+	    call fseti (STDOUT, F_FLUSHNL, YES)
+
+	# Get RFITS parameters.
+	call clgstr ("fits_file", infile, SZ_FNAME)
+	long_header = btoi (clgetb ("long_header"))
+	short_header = btoi (clgetb ("short_header"))
+	len_record = FITS_RECORD
+	old_name = btoi (clgetb ("oldirafname"))
+	make_image = btoi (clgetb ("make_image"))
+
+	# Open the input file list.
+	if (mtfile (infile) == YES) {
+	    inlist = NULL
+	    if (mtneedfileno (infile) == YES)
+	        call clgstr ("file_list", file_list, SZ_LINE)
+	    else
+	        call strcpy ("1", file_list, SZ_LINE)
+	} else {
+	    inlist = fntopnb (infile, NO)
+	    len_inlist = fntlenb (inlist)
+	    if (len_inlist > 0) {
+	        call sprintf  (file_list, SZ_LINE, "1-%d")
+		    call pargi (len_inlist)
+	    } else
+	        call sprintf  (file_list, SZ_LINE, "0")
+	}
+
+	# Decode the ranges string.
+	if (decode_ranges (file_list, range, MAX_RANGES, len_inlist) == ERR)
+	    call error (1, "T_RFITS: Illegal file number list")
+
+	# Open the output file list.
+	if (make_image == YES) {
+	    call clgstr ("iraf_file", outfile, SZ_FNAME)
+	    if (outfile[1] == EOS) {
+		if (old_name == YES)
+		    call mktemp ("tmp", outfile, SZ_FNAME) 
+		else
+		    call error (0, "T_RFITS: Undefined output file name")
+	    }
+	    outlist = fntopnb (outfile, NO)
+	    len_outlist = fntlenb (outlist)
+	    data_type = rft_get_image_type (clgetc ("datatype"))
+	    scale = btoi (clgetb ("scale"))
+	    blank = clgetr ("blank")
+	    offset = clgeti ("offset")
+	} else {
+	    outfile[1] = EOS
+	    outlist = NULL
+	    len_outlist = 1
+	}
+	if ((len_outlist > 1) && (len_outlist != len_inlist))
+	    call error (0,
+	        "T_RFITS: Output and input lists have different lengths")
+
+	# Get the scan size parameter.
+	fe = rft_fe (infile)
+
+	# Read successive FITS files, convert and write into a numbered
+	# succession of output IRAF files.
+
+	file_number = 0
+	while (get_next_number (range, file_number) != EOF) {
+
+	    # Get the input file name.
+	    if (inlist != NULL) {
+		if (fntgfnb (inlist, in_fname, SZ_FNAME) == EOF)
+		    call error (0, "T_RFITS: Error reading input file name")
+	    } else {
+		if (mtneedfileno (infile) == YES)
+		    call mtfname (infile, file_number, in_fname, SZ_FNAME)
+		else
+	            call strcpy (infile, in_fname, SZ_FNAME)
+	    }
+
+	    # Get the output file name.
+	    if (outlist == NULL) {
+		out_fname[1] = EOS
+	    } else if (len_inlist > len_outlist) {
+		if (fntrfnb (outlist, 1, out_fname, SZ_FNAME) == EOF)
+	            call strcpy (outfile, out_fname, SZ_FNAME)
+	        if (len_inlist > 1) {
+		    call sprintf (out_fname[strlen(out_fname)+1], SZ_FNAME,
+		        "%04d")
+		        call pargi (file_number + offset)
+	        }
+	    } else if (fntgfnb (outlist, out_fname, SZ_FNAME) == EOF)
+		call error (0, "T_RFITS: Error reading output file name")
+
+	    # Convert FITS file to the output IRAF file. If EOT is reached
+	    # then exit. If an error is detected then print a warning and
+	    # continue with the next file.
+
+	    iferr (stat = rft_read_fitz (in_fname, out_fname))
+		call erract (EA_FATAL)
+	    if (stat == EOF)
+		break
+	}
+
+	if (inlist != NULL)
+	    call fntclsb (inlist) 
+	if (outlist != NULL)
+	    call fntclsb (outlist) 
+end
+
+
+# RFT_GET_IMAGE_TYPE -- Convert a character to and IRAF image type.
+
+int procedure rft_get_image_type (c)
+
+char	c
+
+int	type_codes[NTYPES], i
+string	types "usilrdx"	
+int	stridx()
+data	type_codes /TY_USHORT, TY_SHORT, TY_INT, TY_LONG, TY_REAL,
+		    TY_DOUBLE, TY_COMPLEX/
+begin
+	i = stridx (c, types)
+	if (i == 0)
+	    return (ERR)
+	else
+	    return (type_codes[stridx(c,types)])
+end
+
+
+# RFT_FE -- Fetch the maximum file size in MB for tape scanning mode.
+
+real procedure rft_fe (file)
+
+char	file[ARB]		# the input file name
+
+pointer	gty
+real	fe
+int	mtfile(), gtygeti()
+pointer	mtcap()
+errchk	gtygeti()
+
+begin
+	if (mtfile (file) == NO)
+	    return (0.0)
+	iferr (gty = mtcap (file))
+	    return (0.0)
+	iferr (fe = gtygeti (gty, "fe"))
+	    fe = 0.0
+	call gtyclose (gty)
+	return (fe)
+end
diff --git a/pkg/obsolete/fits/t_wfits.x b/pkg/obsolete/fits/t_wfits.x
new file mode 100644
index 00000000..13cdd81f
--- /dev/null
+++ b/pkg/obsolete/fits/t_wfits.x
@@ -0,0 +1,216 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include <mach.h>
+include <error.h>
+include <fset.h>
+include "wfits.h"
+
+# T_WFITS -- This procedure converts a series of IRAF image files to
+# FITS image files.
+
+procedure t_wfits ()
+
+char	iraf_files[SZ_FNAME]	# list of IRAF images
+char	fits_files[SZ_FNAME]	# list of FITS files
+bool	newtape			# new or used tape ?
+char	in_fname[SZ_FNAME]	# input file name
+char	out_fname[SZ_FNAME]	# output file name
+
+int	imlist, flist, nimages, nfiles, file_number
+bool	clgetb()
+double	clgetd()
+int	imtopen(), imtlen (), wft_get_bitpix(), clgeti(), imtgetim()
+int	mtfile(), btoi(), fstati(), fntlenb(), fntgfnb(), mtneedfileno()
+int	wft_blkfac(), fntrfnb(), strlen()
+pointer	fntopnb()
+
+include "wfits.com"
+
+begin
+	# Flush on a newline if STDOUT has not been redirected.
+	if (fstati (STDOUT, F_REDIR) == NO)
+	    call fseti (STDOUT, F_FLUSHNL, YES)
+
+	# Open iraf_files template and determine number of files in list.
+	call clgstr ("iraf_files", iraf_files, SZ_FNAME)
+	imlist = imtopen (iraf_files)
+	nimages = imtlen (imlist)
+
+	# Get the wfits parameters.
+	long_header = btoi (clgetb ("long_header"))
+	short_header = btoi (clgetb ("short_header"))
+	make_image = btoi (clgetb ("make_image"))
+
+	# Get the FITS bits per pixel and the FITS logical record size.
+	bitpix = wft_get_bitpix (clgeti ("bitpix"))
+	len_record = FITS_RECORD
+
+	# Get the scaling parameters.
+	scale = btoi (clgetb ("scale"))
+	if (scale == YES) {
+	    if (clgetb ("autoscale"))
+		autoscale = YES
+	    else {
+	        bscale = clgetd ("bscale")
+	        bzero = clgetd ("bzero")
+		autoscale = NO
+	    }
+	} else {
+	    autoscale = NO
+	    bscale = 1.0d0
+	    bzero = 0.0d0
+	}
+
+	# Get the output file name and type (tape or disk). If no tape file
+	# number is given for output, the user is asked if the tape is blank
+	# or contains data. If the tape is blank output begins at BOT,
+	# otherwise at EOT.
+
+	if (make_image == YES) {
+	    call clgstr ("fits_files", fits_files, SZ_FNAME)
+	    if (mtfile (fits_files) == YES) {
+		flist = NULL
+		if (mtneedfileno (fits_files) == YES) {
+	            newtape = clgetb ("newtape")
+		    if (newtape)
+			call mtfname (fits_files, 1, out_fname, SZ_FNAME)
+		    else
+			call mtfname (fits_files, EOT, out_fname, SZ_FNAME)
+	        } else {
+		    call strcpy (fits_files, out_fname, SZ_FNAME)
+		    newtape = false
+		}
+	    } else {
+		flist = fntopnb (fits_files, NO)
+		nfiles = fntlenb (flist)
+		if ((nfiles > 1) && (nfiles != nimages))
+		    call error (0,
+		    "T_WFITS: Input and output lists are not the same length")
+	    }
+	} else {
+	    fits_files[1] = EOS
+	    flist = NULL
+	}
+
+	# Get the fits file blocking factor.
+	blkfac = wft_blkfac (fits_files, clgeti ("blocking_factor"))
+
+	# Loop through the list of input images files.
+
+	file_number = 1
+	while (imtgetim (imlist, in_fname, SZ_FNAME) != EOF) {
+
+	    # Print the id string.
+	    if (long_header == YES || short_header == YES) {
+		call printf ("File %d: %s")
+		    call pargi (file_number)
+		    call pargstr (in_fname)
+	    }
+
+	    # Get the output file name. If single file output to disk, use
+	    # name fits_file. If multiple file output to disk, the file number
+	    # is added to the output file name, if no output name list is
+	    # supplied. If an output name list is supplied then the names
+	    # are extracted one by one from that list.
+
+	    if (make_image == YES) {
+	        if (mtfile (fits_files) == YES) {
+		    if (file_number == 2)
+			call mtfname (out_fname, EOT, out_fname, SZ_FNAME)
+	        } else if (nfiles > 1) {
+		    if (fntgfnb (flist, out_fname, SZ_FNAME) == EOF)
+			 call error (0, "Error reading output file name")
+		} else {
+		    if (fntrfnb (flist, 1, out_fname, SZ_FNAME) == EOF)
+                        call strcpy (fits_files, out_fname, SZ_FNAME)
+		    if (nimages > 1) {
+		        call sprintf (out_fname[strlen(out_fname)+1],
+			    SZ_FNAME, "%04d")
+			    call pargi (file_number)
+		    }
+	        }
+	    }
+
+	    # Write each output file.
+	    iferr (call wft_write_fitz (in_fname, out_fname)) {
+		call printf ("Error writing file: %s\n")
+		    call pargstr (out_fname)
+		call erract (EA_WARN)
+		break
+	    } else
+		file_number = file_number + 1
+	}
+
+	# Close up the input and output lists.
+	call clpcls (imlist)
+	if (flist != NULL)
+	    call fntclsb (flist)
+end
+
+
+# WFT_GET_BITPIX -- This procedure fetches the user determined bitpix or ERR if
+# the bitpix is not one of the permitted FITS types.
+
+int procedure wft_get_bitpix (bitpix)
+
+int	bitpix
+
+begin
+	switch (bitpix) {
+	case FITS_BYTE, FITS_SHORT, FITS_LONG, FITS_REAL, FITS_DOUBLE:
+	    return (bitpix)
+	default:
+	    return (ERR)
+	}
+end
+
+
+# WFT_BLKFAC -- Get the fits tape blocking factor.
+
+int procedure wft_blkfac (file, ublkfac)
+
+char	file[ARB]		# the input file name
+int	ublkfac			# the user supplied blocking factor
+
+int	bs, fb, blkfac
+pointer	gty
+int	mtfile(), mtcap(), gtygeti()
+errchk	mtcap(), gtygeti()
+
+begin
+	# Return a blocking factor of 1 if the file is a disk file.
+	if (mtfile (file) == NO)
+	    return (0)
+
+	# Open the tapecap device entry for the given device, and get
+	# the device block size and default FITS blocking factor
+	# parameters.
+
+	iferr (gty = mtcap (file))
+	    return (max (ublkfac,1))
+	iferr (bs = gtygeti (gty, "bs")) {
+	    call gtyclose (gty)
+	    return (max (ublkfac,1))
+	}
+	iferr (fb = max (gtygeti (gty, "fb"), 1))
+	    fb = 1
+
+	# Determine whether the device is a fixed or variable blocked
+	# device. Set the fits blocking factor to the value of the fb
+	# parameter if device is fixed block or if the user has
+	# requested the default blocking factor. Set the blocking factor
+	# to the user requested value if the device supports variable
+	# blocking factors.
+
+	if (bs == 0) {
+	    if (ublkfac <= 0)
+	        blkfac = fb
+	    else
+		blkfac = ublkfac
+	} else
+	    blkfac = fb
+
+	call gtyclose (gty)
+
+	return (blkfac)
+end
diff --git a/pkg/obsolete/fits/wfits.com b/pkg/obsolete/fits/wfits.com
new file mode 100644
index 00000000..4d2a58f8
--- /dev/null
+++ b/pkg/obsolete/fits/wfits.com
@@ -0,0 +1,15 @@
+# FITS common block
+
+double	bscale		# FITS scaling factor
+double	bzero		# FITS offset factor
+int	bitpix		# Output bits per pixel
+int	len_record	# Record length in FITS bytes
+int	long_header	# Print long header?
+int	short_header	# Print short header?
+int	make_image	# Make a FITS image?
+int	scale		# Scale the data with bzero and bscale?
+int	autoscale	# Allow program to calculate bscale and bzero?
+int	blkfac		# FITS tape blocking factor
+
+common /wfitscom/ bscale, bzero, bitpix, len_record, long_header, short_header,
+		  make_image, scale, autoscale, blkfac
diff --git a/pkg/obsolete/fits/wfits.h b/pkg/obsolete/fits/wfits.h
new file mode 100644
index 00000000..d2a67c34
--- /dev/null
+++ b/pkg/obsolete/fits/wfits.h
@@ -0,0 +1,113 @@
+# WFITS header file
+
+# The basic FITS data structure
+
+define	LEN_FITS	(44 + SZ_FNAME + 1)
+
+define	BSCALE		Memd[P2D($1)]	  # FITS bscale value
+define	BZERO		Memd[P2D($1+2)]	  # FITS bzero value
+define	TAPEMAX		Memd[P2D($1+4)]	  # IRAF tape max
+define	TAPEMIN		Memd[P2D($1+6)]	  # IRAF tape min
+define	IRAFMAX		Memr[P2R($1+8)]	  # IRAF image maximum
+define	IRAFMIN		Memr[P2R($1+9)]	  # IRAF image minimum
+define	BLANK		Meml[P2L($1+10)]  # FITS blank value
+define	FITS_BITPIX	Memi[$1+11]	  # FITS bits per pixel
+define	DATA_BITPIX	Memi[$1+12]	  # Data bits per pixel
+define	SCALE		Memi[$1+13]	  # Scale data?
+define	BLANK_STRING	Memc[P2C($1+19)]  # String containing FITS blank value
+define	TYPE_STRING	Memc[P2C($1+31)]  # String containing IRAF type
+define	IRAFNAME	Memc[P2C($1+41)]  # IRAF file name
+
+
+# Define the FITS record size
+
+define	FITS_RECORD	2880		# Size of standard FITS record (bytes)
+
+# Define the FITS data types
+
+define	FITS_BYTE	8		# Number of bits in a FITS byte
+define	FITS_SHORT	16		# Number of bits in a FITS short
+define	FITS_LONG	32		# Number of bits in a FITS long
+define	FITS_REAL	-32		# Number of bits in a FITS real * -1
+define	FITS_DOUBLE	-64		# Number of bits in a FITS double * -1
+
+# Define the FITS precision in decimal digits
+
+define	BYTE_PREC	3		# Precision of FITS byte
+define	SHORT_PREC	5		# Precision of FITS short
+define	LONG_PREC	10		# Precision of FITS long
+
+# Define the FITS blank data values
+
+define	BYTE_BLANK	0.0d0			# Blank value for a FITS byte
+define	SHORT_BLANK	-3.2768d4		# Blank value for a FITS short
+define	LONG_BLANK	-2.147483648d9		# Blank value for a FITS long
+
+# Define the FITS integer max and min values
+
+define	BYTE_MAX		  2.55d2	# Max value for a FITS byte
+define	BYTE_MIN		   1.0d0	# Min value for a FITS byte
+define	SHORT_MAX	        3.2767d4	# Max value for a FITS short
+define	SHORT_MIN              -3.2767d4	# Min value for a FITS short
+define	LONG_MAX	   2.147483647d9	# Max value for a FITS long
+define	LONG_MIN	  -2.147483647d9	# Min value for a FITS long
+
+# Define the FITS card image parameters
+
+define	LEN_CARD		80	# Length of FITS header card
+define	LEN_KEYWORD		8	# Length of FITS keyword
+define	LEN_NAXIS_KYWRD		5       # Length of the NAXIS keyword string
+define	COL_VALUE		11	# First column of value field
+
+# Mapping of FITS task keywords to IRAF image header keywords
+
+define	NAXIS		IM_NDIM($1)		# Number of dimensions
+define	NAXISN		IM_LEN($1, $2)		# Length of each dimension
+define	OBJECT		IM_TITLE($1)		# Image title
+define	HISTORY		IM_HISTORY($1)		# History
+define	UNKNOWN	Memc[($1+IMU-1)*SZ_STRUCT+1]	# IRAF user area
+
+define	PIXTYPE		IM_PIXTYPE($1)		# Image pixel type
+define	NBPIX		IM_NBPIX($1)		# Number of bad pixels
+define	LIMTIME		IM_LIMTIME($1)		# Last modify limits time
+define	MTIME		IM_MTIME($1)		# Last modify time
+define	CTIME		IM_CTIME($1)		# Create time
+
+define	LEN_USERAREA	28800			# Default user area size
+
+# Set up a structure for the WFITS parameters
+
+# Define the required keywords
+
+define	FIRST_CARD	1		# FITS simple parameter
+define	SECOND_CARD	2		# FITS bitpix parameter
+define	THIRD_CARD	3		# FITS naxis parameter
+
+# Define the optional FITS KEYWORD parameters
+
+define	NOPTIONS	12		# Number of optional keywords
+
+define	KEY_BSCALE	1		# FITS bscale parameter
+define	KEY_BZERO	2		# FITS bzero parameter
+define	KEY_BUNIT	3		# FITS physical units
+define	KEY_BLANK	4		# FITS value of blank pixel
+define	KEY_OBJECT	5		# FITS title string
+define	KEY_ORIGIN	6		# Origin of FITS tape
+define	KEY_DATE	7		# Date the tape was written
+define	KEY_IRAFNAME	8		# Root name of IRAF image
+define	KEY_IRAFMAX	9		# Maximum value of IRAF image
+define	KEY_IRAFMIN	10		# Minimum value of IRAF image
+define	KEY_IRAFBP	11		# Bits per pixel in IRAF image
+define	KEY_IRAFTYPE	12		# IRAF image data type
+
+define	LEN_STRING	8		# Minimum length of a string parameter
+define	LEN_OBJECT	63		# Maximum length of string parameter
+define	LEN_ALIGN	18		# Maximum length for aligning parameter
+define	LEN_ORIGIN	9		# Length of origin string
+define	LEN_BLANK	11		# Length of the blank string
+define	NDEC_REAL	7		# Precision of real data
+define	NDEC_DOUBLE	11		# Precision of double precision data
+
+# Miscellaneous
+
+define	CENTURY		1900