Repatch (from linux) of OSX IRAF

1 files changed, 1202 insertions, 0 deletions

diff --git a/pkg/plot/t_implot.x b/pkg/plot/t_implot.x
new file mode 100644
index 00000000..fbf8f00a
--- /dev/null
+++ b/pkg/plot/t_implot.x
@@ -0,0 +1,1202 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<error.h>
+include	<ctype.h>
+include	<imhdr.h>
+include	<mach.h>
+include	<gset.h>
+include	<mwset.h>
+
+define	SEGLEN		10
+define	SZ_PLOTTITLE	512
+define	KEYSFILE	"lib$scr/implot.key"
+define	MAX_COLORS	8
+
+
+# IMPLOT -- Image plot program.  An interactive, cursor driven program for
+# plotting lines and columns of images.  This is an early version of the program
+# lacking averaging and interaction with the image display cursor.
+#
+#	implot (image [,line])
+#
+# Keystrokes:
+#
+#	?		help
+#	a		plot the average of a range of lines of columns
+#	c		plot column at position of cursor
+#	e		expand plot by marking corners of new window
+#	j		move down
+#	k		move up
+#	l		plot line at position of cursor
+#	m		previous image
+#	n		next image
+#	p		measure profile (mark region and baseline with 2 pos)
+#	o		overplot next vector
+#	r		redraw
+#	s		print statistics on a region
+#	/		scroll status line
+#
+#
+# In addition to the above keystrokes, the following ':' escapes are recognized
+# by the program:
+#
+#	:a N		set number of lines or columns to average
+#	:c M [N]	plot column M or avg of M to N
+#	:f format	set label format
+#	:i imagename	open a new image for input
+#	:l M [N]	plot line M or avg of M to N
+#	:o		overplot
+#	:log+,log-	enable, disable log scaling in Y
+#	:step N		set step size for j,k
+#	:solid		overplot with solid, not dashed, lines
+#	:mono		disable coloring of overplotted vectors
+#	:x x1 x2	fix plot window in X (no args to unfix)
+#	:y y1 y2	fix plot window in Y (no args to unfix)
+#	:w wcstype	change world coordinate type
+
+procedure t_implot()
+
+int	list
+char	image[SZ_FNAME]
+char	wcstype[SZ_FNAME]
+char	xlabel[SZ_FNAME]
+char	format[SZ_FNAME]
+char	fmt[SZ_FNAME]
+char	command[SZ_FNAME]
+char	device[SZ_FNAME]
+
+int	xnticks, ynticks
+bool	overplot, lineplot, logscale, erase, rescale[2], p_rescale[2]
+int	key, wcs, ip, i1, i2, n, linetype, color, nltypes, linestep, navg
+int	npix, nlines, ncols, line, col, shift, step, p_navg, sline, nim, index
+real	x, y, px, py, qx, qy, x1, x2, y1, y2
+real	median, mean, sigma, sum
+pointer	im, mw, ct, gp, xold, yold, xnew, ynew, sl, ptr
+
+real	asumr(), amedr(), plt_iformatr()
+int	clgeti(), clgcur(), ctoi(), ctor(), ggeti(), imtlen(), imaccess()
+pointer	gopen(), immap(), mw_openim(), mw_sctran(), sl_getstr()
+pointer	imtopenp(), imtrgetim()
+errchk	mw_sctran
+
+define	line_ 91
+define	col_ 92
+define	replotline_ 93
+define	replotcol_ 94
+define	nextim_ 95
+define	quit_ 96
+string	bell "\007"
+string	again "again:"
+
+begin
+	list = imtopenp ("image")
+	call clgstr ("device", device, SZ_FNAME)
+	gp = gopen (device, NEW_FILE, STDGRAPH)
+	call clgstr ("wcs", wcstype, SZ_FNAME)
+
+	if (clgeti ("$nargs") > 1)
+	    line = clgeti ("line")
+	else
+	    line = INDEFI
+
+	p_rescale[1] = true
+	rescale[1]   = true
+	p_rescale[2] = true
+	rescale[2]   = true
+
+	logscale  = false
+	overplot  = false
+	lineplot  = true
+	erase     = false
+	xnticks	  = 5
+	ynticks	  = 5
+	format[1] = EOS
+
+	linestep  = 1
+	linetype  = 1
+	color     = 1
+	nltypes   = ggeti (gp, "lt")
+	p_navg	  = 1
+	navg	  = 1
+	step	  = clgeti ("step")
+
+	# Loop through the images.  Currently this loop is not actually
+	# used and instead the 'm' and 'n' keys explicitly change the
+	# image.  The 'q' key exits the loop regardless of the position
+	# of the list.
+
+	nim = imtlen (list)
+	index = 1
+	while (imtrgetim (list, index, image, SZ_FNAME) != EOF) {
+	    iferr {
+		im = NULL; mw = NULL; sl = NULL
+
+		if (imaccess (image, READ_ONLY) == YES) {
+		    ptr = immap (image, READ_ONLY, 0); 
+		    im = ptr
+	   
+		} else {
+		    call eprintf ("Error opening image '%s'")
+			call pargstr (image)
+		    goto nextim_
+		}
+
+		ptr = mw_openim (im); mw = ptr
+		call mw_seti (mw, MW_USEAXMAP, NO)
+
+		ct = mw_sctran (mw, "logical", wcstype, 0)
+
+		ncols = IM_LEN(im,1)
+		if (IM_NDIM(im) <= 0)
+		    call error (1, "image has no pixels")
+		else if (IM_NDIM(im) > 1)
+		    nlines = IM_LEN(im,2)
+		else
+		    nlines = 1
+
+		if (IS_INDEFI(line))
+		    line = max(1, min(nlines, (nlines + 1) / 2))
+
+		if (IS_INDEFI(step) || step < 1)
+		    step = max (1, nlines / 10)
+
+		npix = max (ncols, nlines)
+		call malloc (xold, npix, TY_REAL)
+		call malloc (yold, npix, TY_REAL)
+		call malloc (xnew, npix, TY_REAL)
+		call malloc (ynew, npix, TY_REAL)
+
+		if (!overplot)
+		    call gclear (gp)
+		call imp_getvector (im, mw, ct, wcstype, xnew, ynew, xlabel,
+		    fmt, line, navg, lineplot)
+		if (format[1] == '%')
+		    call strcpy (format, fmt, SZ_FNAME)
+		npix = ncols
+
+		call gsets (gp, G_XTICKFORMAT, fmt)
+		if (xnticks >= 0)
+		    call gseti (gp, G_XNMAJOR, xnticks)
+		if (ynticks >= 0)
+		    call gseti (gp, G_YNMAJOR, ynticks)
+		call gseti (gp, G_NMINOR, 0)
+		call gseti (gp, G_PLTYPE, 1)
+		call gseti (gp, G_PLCOLOR, 1)
+
+		call imp_plotvector (gp, im, Memr[xnew], Memr[ynew], ncols,
+		    nlines, real(line), navg, lineplot, rescale, image, xlabel)
+		overplot = false
+
+		call sl_init (sl, 1)
+		while (clgcur ("coords", x, y, wcs, key, command, SZ_FNAME) !=
+		    EOF) {
+		    if (key == 'q')
+quit_			break
+
+		    switch (key) {
+		    case 'a':
+			# Plot the average over a range of lines or columns
+			# marked interactively with the cursor.
+
+			x1 = x; y1 = y
+			call printf (again)
+			if (clgcur ("gcur", x2, y2, wcs, key, command,
+			    SZ_FNAME) == EOF)
+			    next
+
+			if (abs(x2-x1) > abs(y2-y1)) {
+			    # Range is in X.
+
+			    navg = abs (x2 - x1) + 1
+			    if (lineplot) {
+				col = nint (min (x1, x2))
+				goto col_
+			    } else {
+				line = nint (min (x1, x2))
+				goto line_
+			    }
+
+			} else {
+			    # Range is in Y.
+
+			    if (lineplot) {
+				call imp_tran (gp, x1, y1, x1, y1, Memr[xnew],
+				    ncols, nlines)
+				call imp_tran (gp, x2, y2, x2, y2, Memr[xnew],
+				    ncols, nlines)
+				navg = abs (y2 - y1) + 1
+				line = nint (min (y1, y2))
+				goto line_
+			    } else {
+				call imp_tran (gp, x1, y1, x1, y1, Memr[xnew],
+				    nlines, ncols)
+				call imp_tran (gp, x2, y2, x2, y2, Memr[xnew],
+				    nlines, ncols)
+				navg = abs (y2 - y1) + 1
+				col = nint (min (y1, y2))
+				goto col_
+			    }
+			}
+
+		    case 'j', 'k':
+			# Move viewport into image up (k) or down (j).  This
+			# is done by erasing the old data vector and drawing
+			# a new one.
+
+			erase = true
+			navg = p_navg
+			overplot = true
+			call amovr (Memr[xnew], Memr[xold], npix)
+			call amovr (Memr[ynew], Memr[yold], npix)
+
+			shift = step
+			if (key == 'j')
+			    shift = -shift
+			    
+			if (lineplot) {
+			    line = line + shift
+			    goto line_
+			} else {
+			    col  = col  + shift
+			    goto col_
+			}
+
+		    case 'l':
+			# Plot a line.
+			if (lineplot) {
+			    call imp_tran (gp, x, y, px, py, Memr[xnew], ncols,
+				nlines)
+			    line = max(1, min(nlines, nint(py)))
+			} else {
+			    call imp_tran (gp, x, y, px, py, Memr[xnew], nlines,
+				ncols)
+			    line = max(1, min(nlines, nint(px)))
+			}
+			navg = p_navg
+			line = line - (navg - 1) / 2
+line_
+			lineplot = true
+			line = max(1, min(nlines, line))
+			call imp_getvector (im, mw, ct, wcstype, xnew, ynew,
+			    xlabel, fmt, line, navg, lineplot)
+			if (format[1] == '%')
+			    call strcpy (format, fmt, SZ_FNAME)
+			npix = ncols
+replotline_
+			if (overplot) {
+			    if (erase) {
+				# Erase old vector and replace it with new
+				# vector.
+
+				call imp_redraw (gp, Memr[xold], Memr[yold],
+				    Memr[xnew], Memr[ynew], npix)
+				erase = false
+
+			    } else {
+				# Overplot new vector. 
+
+				linetype = linetype + linestep
+				if (linetype > nltypes)
+				    linetype = 1
+				call gseti (gp, G_PLTYPE, linetype)
+
+				color = color + 1
+				if (color > MAX_COLORS)
+				    color = 1
+				call gseti (gp, G_PLCOLOR, color)
+
+				call gpline (gp, Memr[xnew], Memr[ynew], ncols)
+			    }
+
+			    call imp_markpos (gp, line, nlines)
+			    overplot = false
+
+			} else {
+			    call gclear (gp)
+			    call gsets (gp, G_XTICKFORMAT, fmt)
+			    if (logscale)
+				call gseti (gp, G_YTRAN, GW_LOG)
+			    call gseti (gp, G_NMINOR, 0)
+			    if (xnticks >= 0)
+				call gseti (gp, G_XNMAJOR, xnticks)
+			    if (ynticks >= 0)
+				call gseti (gp, G_YNMAJOR, ynticks)
+			    linetype = 1
+			    color = 1
+			    call imp_plotvector (gp, im, Memr[xnew], Memr[ynew],
+				ncols, nlines, real(line), navg, lineplot,
+				rescale, image, xlabel)
+			    rescale[1] = p_rescale[1]
+			    rescale[2] = p_rescale[2]
+			}
+
+		    case 'm', 'n':
+			if (key == 'm') {
+			    if (index > 1)
+				index = index - 1
+			    else
+				next
+			} else if (key == 'n') {
+nextim_		    	if (index < nim)
+				index = index + 1
+			    else
+				next
+			}
+
+			if (imtrgetim (list, index, command, SZ_FNAME) == EOF)
+			    break
+
+			# Open a different image.
+			call mw_close (mw)
+			call imunmap (im)
+
+			iferr (im = immap (command, READ_ONLY, 0)) {
+			    call erract (EA_WARN)
+			    im = immap (image, READ_ONLY, 0)
+			    mw = mw_openim (im)
+			    call mw_seti (mw, MW_USEAXMAP, NO)
+			    ct = mw_sctran (mw, "logical", wcstype, 0)
+			    next
+			}
+
+			if (IM_NDIM(im) <= 0) {
+			    call eprintf ("image has no pixels\n")
+			    im = immap (image, READ_ONLY, 0)
+			    mw = mw_openim (im)
+			    call mw_seti (mw, MW_USEAXMAP, NO)
+			    ct = mw_sctran (mw, "logical", wcstype, 0)
+			    next
+
+			} else {
+			    mw = mw_openim (im)
+			    call mw_seti (mw, MW_USEAXMAP, NO)
+			    ct = mw_sctran (mw, "logical", wcstype, 0)
+
+			    ncols = IM_LEN(im,1)
+			    if (IM_NDIM(im) > 1)
+				nlines = IM_LEN(im,2)
+			    else {
+				lineplot = true
+				nlines = 1
+			    }
+
+			    npix   = max (ncols, nlines)
+			    call strcpy (command, image, SZ_FNAME)
+			    call realloc (xold, npix, TY_REAL)
+			    call realloc (yold, npix, TY_REAL)
+			    call realloc (xnew, npix, TY_REAL)
+			    call realloc (ynew, npix, TY_REAL)
+
+			    if (lineplot)
+				goto line_
+			    else
+				goto col_
+			}
+
+		    case 'c':
+			# Plot a column.
+			if (lineplot) {
+			    call imp_tran (gp, x, y, px, py, Memr[xnew], ncols,
+				nlines)
+			    col = max(1, min(ncols, nint(px)))
+			} else {
+			    call imp_tran (gp, x, y, px, py, Memr[xnew], nlines,
+				ncols)
+			    col = max(1, min(ncols, nint(py)))
+			}
+			navg = p_navg
+			col = col - (navg - 1) / 2
+col_
+			if (nlines == 1) {
+			    call printf (bell)
+			    next
+			}
+			lineplot = false
+			col = max(1, min(ncols, col))
+			call imp_getvector (im, mw, ct, wcstype, xnew, ynew,
+			    xlabel, fmt, col, navg, lineplot)
+			if (format[1] == '%')
+			    call strcpy (format, fmt, SZ_FNAME)
+			npix = nlines
+replotcol_
+			if (overplot) {
+			    if (erase) {
+				# Erase old vector and replace it with new
+				# vector.
+
+				call imp_redraw (gp, Memr[xold], Memr[yold],
+				    Memr[xnew], Memr[ynew], npix)
+				erase = false
+
+			    } else {
+				linetype = linetype + linestep
+				if (linetype > nltypes)
+				    linetype = 1
+				call gseti (gp, G_PLTYPE, linetype)
+
+				color = color + 1
+				if (color > MAX_COLORS)
+				    color = 1
+				call gseti (gp, G_PLCOLOR, color)
+
+				call gpline (gp, Memr[xnew], Memr[ynew], nlines)
+			    }
+
+			    call imp_markpos (gp, col, ncols)
+			    overplot = false
+
+			} else {
+			    call gclear (gp)
+			    call gsets (gp, G_XTICKFORMAT, fmt)
+			    if (logscale)
+				call gseti (gp, G_YTRAN, GW_LOG)
+			    call gseti (gp, G_NMINOR, 0)
+			    if (xnticks >= 0)
+				call gseti (gp, G_XNMAJOR, xnticks)
+			    if (ynticks >= 0)
+				call gseti (gp, G_YNMAJOR, ynticks)
+			    linetype = 1
+			    color = 1
+			    call imp_plotvector (gp, im, Memr[xnew], Memr[ynew],
+				nlines, ncols, real(col), navg, lineplot,
+				rescale, image, xlabel)
+			    rescale[1] = p_rescale[1]
+			    rescale[2] = p_rescale[2]
+			}
+
+		    case 'e':
+			# Expand plot by marking corners of new window.  We are
+			# called with the coords of the lower left corner.
+
+			x1 = x; y1 = y
+			call printf (again)
+			if (clgcur ("gcur", x2, y2, wcs, key, command,
+			    SZ_FNAME) == EOF)
+			    next
+
+			rescale[1] = false
+			rescale[2] = false
+			p_rescale[1] = true
+			p_rescale[2] = true
+
+			# If the cursor moved only in X, with negligible range
+			# in Y, expand only in X.  Do the comparisons in NDC
+			# space to avoid scaling problems.
+
+			call gctran (gp, x1, y1, px, py, wcs, 0)
+			call gctran (gp, x2, y2, qx, qy, wcs, 0)
+
+			if (abs (py - qy) < .01) {
+			    y1 = INDEF;  y2 = INDEF
+			    rescale[2] = true
+			}
+			call imp_swind (x1, x2, y1, y2)
+
+			if (lineplot)
+			    goto replotline_
+			else
+			    goto replotcol_
+
+		    case 'o':
+			overplot = true
+
+		    case 'r':
+			if (lineplot)
+			    goto replotline_
+			else
+			    goto replotcol_
+
+		    case 'p':
+			# Profile analysis.
+			x1 = x
+			y1 = y
+			call printf (again)
+			if (clgcur ("gcur", x2, y2, wcs, key, command,
+			    SZ_FNAME) == EOF)
+			    next
+
+			call imp_profile (gp, Memr[xnew], Memr[ynew], npix,
+			    x1, y1, x2, y2, sl, sline)
+			call printf (Memc[sl_getstr(sl,sline)])
+
+		    case 's':
+			# Statistics.
+			x1 = x
+			call printf (again)
+			if (clgcur ("gcur", x2, y, wcs, key, command,
+			    SZ_FNAME) == EOF)
+			    next
+
+			i1 = max(1, min(npix, nint(x1)))
+			i2 = max(1, min(npix, nint(x2)))
+			if (i1 > i2) {
+			    n = i1
+			    i1 = i2
+			    i2 = n
+			} else if (i1 == i2)
+			    i2 = i1 + 1
+
+			n = i2 - i1 + 1
+			call aavgr (Memr[ynew+i1-1], n, mean, sigma)
+			median = amedr (Memr[ynew+i1-1], n)
+			sum = asumr (Memr[ynew+i1-1], n)
+
+			call sl_init (sl, 1)
+			call sprintf (Memc[sl_getstr(sl,1)], SZ_LINE,
+			    "median=%g, mean=%g, rms=%g, sum=%g, npix=%d\n")
+			    call pargr (median)
+			    call pargr (mean)
+			    call pargr (sigma)
+			    call pargr (sum)
+			    call pargi (n)
+			sline = 1
+			call printf (Memc[sl_getstr(sl,sline)])
+
+		    case ' ':
+			# Print cursor coordinates.
+			call sl_init (sl, 1)
+			if (lineplot) {
+			    call imp_tran (gp, x, y, px, py, Memr[xnew], ncols,
+				nlines)
+			    col = px
+			    call plt_wcscoord (im, mw, ct, wcstype, format, col,
+				line, Memr[ynew+col-1], Memc[sl_getstr(sl,1)],
+				SZ_LINE)
+			} else {
+			    call imp_tran (gp, x, y, px, py, Memr[xnew], nlines,
+				ncols)
+			    line = px
+			    call plt_wcscoord (im, mw, ct, wcstype, format, col,
+				line, Memr[ynew+line-1], Memc[sl_getstr(sl,1)],
+				SZ_LINE)
+			}
+			sline = 1
+			call printf (Memc[sl_getstr(sl,sline)])
+
+		    case '?':
+			# Print command summary.
+			call gpagefile (gp, KEYSFILE, "implot cursor commands")
+
+		    case ':':
+			# Command mode.
+			for (ip=1;  IS_WHITE (command[ip]);  ip=ip+1)
+			    ;
+			if (command[ip] == 'o') {
+			    overplot = true
+			    ip = ip + 1
+			}
+
+			switch (command[ip]) {
+			case 'a':
+			    # Set number of lines or columns to average.
+			    ip = ip + 1
+			    if (ctoi (command, ip, p_navg) <= 0) {
+				call printf (bell)
+				p_navg = 1
+			    }
+
+			case 'i':
+			    # Open a different image.
+			    call mw_close (mw)
+			    call imunmap (im)
+			    ip = ip + 1
+			    while (IS_WHITE (command[ip]))
+				ip = ip + 1
+
+			    iferr (im = immap (command[ip], READ_ONLY, 0)) {
+				call erract (EA_WARN)
+				im = immap (image, READ_ONLY, 0)
+				mw = mw_openim (im)
+				call mw_seti (mw, MW_USEAXMAP, NO)
+				ct = mw_sctran (mw, "logical", wcstype, 0)
+
+			    } else if (IM_NDIM(im) <= 0) {
+				call eprintf ("image has no pixels\n")
+				im = immap (image, READ_ONLY, 0)
+				mw = mw_openim (im)
+				call mw_seti (mw, MW_USEAXMAP, NO)
+				ct = mw_sctran (mw, "logical", wcstype, 0)
+
+			    } else {
+				mw = mw_openim (im)
+				call mw_seti (mw, MW_USEAXMAP, NO)
+				ct = mw_sctran (mw, "logical", wcstype, 0)
+
+				ncols = IM_LEN(im,1)
+				if (IM_NDIM(im) > 1)
+				    nlines = IM_LEN(im,2)
+				else
+				    nlines = 1
+
+				npix   = max (ncols, nlines)
+				call strcpy (command[ip], image, SZ_FNAME)
+				call realloc (xold, npix, TY_REAL)
+				call realloc (yold, npix, TY_REAL)
+				call realloc (xnew, npix, TY_REAL)
+				call realloc (ynew, npix, TY_REAL)
+			    }
+
+			case 'w':
+			    # Change wcs type.
+			    call mw_ctfree (ct)
+			    ip = ip + 1
+			    while (IS_WHITE (command[ip]))
+				ip = ip + 1
+
+			    iferr {
+				ct = mw_sctran (mw, "logical", command[ip], 0)
+				call strcpy (command[ip], wcstype, SZ_FNAME)
+			    } then {
+				call erract (EA_WARN)
+				ct = mw_sctran (mw, "logical", wcstype, 0)
+			    } else {
+				# Only replot if WCS command succeeds,
+				# otherwise the error message is lost.
+				if (lineplot)
+				    goto line_
+				else
+				    goto col_
+			    }
+
+			case 'f':
+			    # Change label format.
+			    ip = ip + 1
+			    while (IS_WHITE (command[ip]))
+				ip = ip + 1
+			    if (command[ip] == '%') {
+				call strcpy (command[ip], format, SZ_FNAME)
+				call strcpy (format, fmt, SZ_FNAME)
+				if (lineplot)
+				    goto replotline_
+				else
+				    goto replotcol_
+			    } else if (format[1] == '%') {
+				call strcpy (command[ip], format, SZ_FNAME)
+				if (lineplot)
+				    goto line_
+				else
+				    goto col_
+			    }
+
+			case 'l':
+			    if (command[ip+1] != 'o') {
+				# Plot a line.
+				ip = ip + 1
+				if (ctoi (command, ip, i1) <= 0) {
+				    call printf (bell)
+				    next
+				} else if (ctoi (command, ip, i2) <= 0) {
+				    line = max(1, min(nlines, i1))
+				    navg = p_navg
+				    line = line - (navg - 1) / 2
+				    goto line_
+				} else {
+				    i1 = max(1, min(nlines, i1))
+				    i2 = max(1, min(nlines, i2))
+				    line = min (i1, i2)
+				    navg = max (1, abs (i2 - i1) + 1)
+				    goto line_
+				}
+			    } else {
+				# Enable/disable log scaling.
+				while (IS_ALPHA(command[ip]))
+				    ip = ip + 1
+				logscale = (command[ip] == '+')
+			    }
+
+			case 'c':
+			    # Plot a column.
+			    ip = ip + 1
+			    if (ctoi (command, ip, i1) <= 0) {
+				call printf (bell)
+				next
+			    } else if (ctoi (command, ip, i2) <= 0) {
+				col = max(1, min(ncols, i1))
+				navg = p_navg
+				col = col - (navg - 1) / 2
+				goto col_
+			    } else {
+				i1 = max(1, min(ncols, i1))
+				i2 = max(1, min(ncols, i2))
+				col  = min (i1, i2)
+				navg = max (1, abs (i2 - i1) + 1)
+				goto col_
+			    }
+
+			case 's':
+			    if (command[ip+1] == 'o') {
+				# Use only linetype=1 (solid).
+				linetype = 1
+				linestep = 0
+				color = 1
+			    } else {
+				# Set step size.
+				while (IS_ALPHA (command[ip]))
+				    ip = ip + 1
+				if (ctoi (command, ip, step) <= 0) {
+				    call printf (bell)
+				    step = 1
+				}
+			    }
+			
+			case 'x':
+			    # Fix window in X and replot vector.  If no args
+			    # are given, unfix the window.
+
+			    ip = ip + 1
+			    if (ctor (command, ip, x1) <= 0) {
+				rescale[1]   = true
+				p_rescale[1] = true
+			    } else if (ctor (command, ip, x2) <= 0) {
+				call printf (bell)
+			    } else {
+				x1 = plt_iformatr (x1, fmt)
+				x2 = plt_iformatr (x2, fmt)
+				call imp_swind (x1, x2, INDEF, INDEF)
+				rescale[1]   = false
+				p_rescale[1] = false
+			    }
+
+			    if (lineplot)
+				goto replotline_
+			    else
+				goto replotcol_
+			
+			case 'y':
+			    # Fix window in Y and replot vector.  If no args
+			    # are given, unfix the window.
+
+			    ip = ip + 1
+			    if (ctor (command, ip, y1) <= 0) {
+				rescale[2]   = true
+				p_rescale[2] = true
+			    } else if (ctor (command, ip, y2) <= 0) {
+				call printf (bell)
+			    } else {
+				y1 = plt_iformatr (y1, fmt)
+				y2 = plt_iformatr (y2, fmt)
+				call imp_swind (INDEF, INDEF, y1, y2)
+				p_rescale[2] = false
+				rescale[2]   = false
+			    }
+
+			    if (lineplot)
+				goto replotline_
+			    else
+				goto replotcol_
+
+			case 'n':
+			    ip = ip + 1
+			    if (command[ip] == 'x') {
+				while (IS_ALPHA(command[ip]))
+				    ip = ip + 1
+				if (ctoi (command, ip, xnticks) <= 0)
+				    xnticks = -1
+			    } else if (command[ip] == 'y') {
+				while (IS_ALPHA(command[ip]))
+				    ip = ip + 1
+				if (ctoi (command, ip, ynticks) <= 0)
+				    ynticks = -1
+			    } else
+				call printf (bell)
+
+			default:
+			    call printf (bell)
+			}
+
+		    case '/':
+			# Scroll or rewrite the status line.
+
+			sline = sline + 1
+			call printf (Memc[sl_getstr(sl,sline)])
+
+		    default:
+			call printf (bell)
+		    }
+		}
+	    } then
+		call erract (EA_WARN)
+
+	    call mfree (xnew, TY_REAL)
+	    call mfree (ynew, TY_REAL)
+	    call mfree (xold, TY_REAL)
+	    call mfree (yold, TY_REAL)
+	    if (sl != NULL)
+		call sl_free (sl)
+
+	    if (mw != NULL)
+		call mw_close (mw)
+	    if (im != NULL)
+		call imunmap (im)
+
+	    if (key == 'q')
+		break
+	}
+
+	call gclose (gp)
+	call imtclose (list)
+end
+
+
+# IMP_GETVECTOR -- Get a data vector, i.e., line or column or average of
+# lines and columns.
+
+procedure imp_getvector (im, mw, ct, wcstype, x, y, xlabel, format, linecol,
+	navg, lineplot)
+
+pointer	im			# image descriptor
+pointer	mw			# mwcs descriptor
+pointer	ct			# coordinate descriptor
+char	wcstype[ARB]		# WCS type
+pointer	x, y			# output vector
+char	xlabel[SZ_FNAME]	# WCS label
+char	format[SZ_FNAME]	# WCS format
+int	linecol			# line or column number
+int	navg			# number of lines or columns to be averaged
+bool	lineplot		# true if line is to be extracted
+
+real	norm
+pointer	sp, axvals, buf, off
+int	x1, x2, y1, y2
+int	nx, ny, width, i, ndim
+real	asumr()
+pointer	imgl2r(), imgs2r(), imgl1r(), imgs1r()
+errchk	imgl2r, imgs2r, imgl1r, imgs1r, plt_wcs
+
+begin
+	call smark (sp)
+	call salloc (axvals, IM_NDIM(im), TY_REAL)
+
+	call strcpy (wcstype, xlabel, SZ_FNAME)
+
+	ndim = IM_NDIM(im)
+	nx = IM_LEN(im,1)
+	if (ndim > 1)
+	    ny = IM_LEN(im,2)
+	else
+	    ny = 1
+	call amovkr (1., Memr[axvals], ndim)
+
+	if (lineplot) {
+	    # Extract a line vector.
+
+	    x1 = 1
+	    x2 = nx
+	    y1 = max(1, min (ny, linecol))
+	    y2 = max(1, min (ny, linecol + navg - 1))
+
+	    Memr[axvals+1] = (y1 + y2) / 2.
+	    call plt_wcs (im, mw, ct, 1, Memr[axvals], real(x1), real(x2),
+		Memr[x], nx, xlabel, format, SZ_FNAME)
+
+	    if (ndim == 1)
+		call amovr (Memr[imgl1r(im)], Memr[y], nx)
+
+	    else {
+		# Compute sum.
+		call aclrr (Memr[y], nx)
+		do i = y1, y2
+		    call aaddr (Memr[imgl2r(im,i)], Memr[y], Memr[y], nx)
+
+		# Normalize.
+		width = y2 - y1 + 1
+		if (width > 1)
+		    call amulkr (Memr[y], 1. / width, Memr[y], nx)
+	    }
+
+	} else {
+	    # Extract a column vector.
+
+	    x1 = max(1, min(nx, linecol))
+	    x2 = max(1, min(nx, linecol + navg - 1))
+	    y1 = 1
+	    y2 = ny
+
+	    Memr[axvals] = (x1 + x2) / 2.
+	    call plt_wcs (im, mw, ct, 2, Memr[axvals], real(y1), real(y2),
+		Memr[x], ny, xlabel, format, SZ_LINE)
+
+	    width = x2 - x1 + 1
+	    norm  = 1.0 / real(width)
+
+	    if (width > 1) {
+		call aclrr (Memr[y], ny)
+		do i = y1, y2 {
+		    if (ndim == 1) {
+			buf = imgs1r (im, x1, x2)
+			off = buf
+		    } else if (nx > 1024) {
+			buf = imgs2r (im, x1, x2, i, i)
+			off = buf
+		    } else {
+			buf = imgl2r (im, i)
+			off = buf + x1 - 1
+		    }
+		    Memr[y+i-1] = asumr (Memr[off], width) * norm
+		}
+	    } else {
+		buf = imgs2r (im, x1, x2, y1, y2)
+		call amovr (Memr[buf], Memr[y], ny)
+	    }
+	}
+
+	call sfree (sp)
+end
+
+
+# IMP_PLOTVECTOR -- Plot a line or column vector.
+
+procedure imp_plotvector (gp, im, xv, yv, nx, ny, y, navg, lineplot, rescale,
+	image, xlabel)
+
+pointer	gp			# graphics descriptor
+pointer	im			# image descriptor
+real	xv[ARB]			# coordinate vector
+real	yv[ARB]			# data  vector
+int	nx			# number of pixels in vector
+int	ny			# number of pixels on plot-Y axis
+real	y			# position on plot Y-axis
+int	navg			# number of lines or columns averaged
+bool	lineplot		# are we plotting a line or a column
+bool	rescale[2]		# rescale plot
+char	image[ARB]		# image name
+char	xlabel[ARB]		# X label
+
+real	junkr
+int	i, i1, i2, npix, maxch
+pointer	sp, ip, plot_title, op
+bool	fp_equalr()
+
+real	x1, x2, y1, y2
+common	/implcom/ x1, x2, y1, y2
+
+begin
+	call smark (sp)
+	call salloc (plot_title, SZ_PLOTTITLE, TY_CHAR)
+
+	# Format the plot title, starting with the system banner.
+	call sysid (Memc[plot_title], SZ_PLOTTITLE)
+	for (op=plot_title;  Memc[op] != '\n' && Memc[op] != EOS;  op=op+1)
+	    ;
+	Memc[op] = '\n'
+	op = op + 1
+	maxch = SZ_PLOTTITLE - (op - plot_title)
+	# Format the remainder of the plot title.
+
+	if (IM_LEN(im,2) <= 1) {
+	    # Plot of a one-dim image.
+	    call strcpy (IM_TITLE(im), Memc[op], maxch)
+
+	} else if (navg > 1) {
+	    call sprintf (Memc[op], maxch, "Average of %s %d to %d of %s\n%s")
+		if (lineplot) {
+		    call pargstr ("lines")
+		    npix = IM_LEN(im,2)
+		} else {
+		    call pargstr ("columns")
+		    npix = IM_LEN(im,1)
+		}
+
+		i1 = max (1, min (npix, nint (y)))
+		i2 = max (1, min (npix, nint (y) + navg - 1))
+		call pargi (i1)
+		call pargi (i2)
+		call pargstr (image)
+		call pargstr (IM_TITLE(im))
+
+	} else {
+	    call sprintf (Memc[op], maxch, "%s %d of %s\n%s")
+		if (lineplot)
+		    call pargstr ("Line")
+		else
+		    call pargstr ("Column")
+		call pargi (nint(y))
+		call pargstr (image)
+		call pargstr (IM_TITLE(im))
+	}
+
+
+	# Delete trailing newline and any whitespace from image title string.
+	# Trailing whitespace causes plot title to not be drawn centered on
+	# plot.
+
+	for (ip=plot_title;  Memc[ip] != EOS;  ip=ip+1)
+	    ;
+	ip = ip - 1
+	if (ip > plot_title && Memc[ip] == '\n')
+	    ip = ip - 1
+	while (ip > plot_title && IS_WHITE(Memc[ip]))
+	    ip = ip - 1
+	Memc[ip+1] = EOS
+
+	# Autoscale the plot in X and or Y if so indicated.
+	if (rescale[1])
+	    call gascale (gp, xv, nx, 1)
+	else
+	    call gswind (gp, x1, x2, INDEF, INDEF)
+
+	call ggwind (gp, x1, x2, junkr, junkr)
+	junkr = min (x1, x2)
+	for (i1=1; i1<nx && xv[i1] < junkr; i1=i1+1)
+	    ;
+	junkr = max (x1, x2)
+	for (i2=nx; i2>1 && xv[i2] > junkr; i2=i2-1)
+	    ;
+	if (i2 < i1) {
+	    i = i1
+	    i1 = i2
+	    i2 = i
+	}
+	npix = max (1, i2 - i1 + 1)
+
+	if (rescale[2]) {
+	    if (npix < 2)
+		call gascale (gp, yv[i1], nx, 2)
+	    else
+		call gascale (gp, yv[i1], npix, 2)
+	} else
+	    call gswind (gp, INDEF, INDEF, y1, y2)
+
+	call ggwind (gp, x1, x2, y1, y2)
+
+	# If the image is two dimensional plot the position within the image
+	# of the plotted vector on the plot-Y axis (which may refer to either
+	# X or Y on the image).
+
+	if (IM_LEN(im,2) > 1) {
+	    # Draw all but right axes.
+	    if (fp_equalr (y1, y2)) {
+		y1 = 0.99 * y1
+		y2 = 1.01 * y2
+		call gswind (gp, INDEF, INDEF, y1, y2)
+	    }
+	    call gseti (gp, G_YDRAWAXES, 1)
+	    call glabax (gp, Memc[plot_title], xlabel, "")
+
+	    # Draw right axis (pixel scale)
+	    call ggwind (gp, x1, x2, y1, y2)
+	    call gswind (gp, 1., real (nx), 1., real (ny))
+	    call gseti (gp, G_XDRAWAXES, 0)
+	    call gseti (gp, G_YDRAWAXES, 2)
+	    call glabax (gp, "", "", "")
+	    call gswind (gp, x1, x2, y1, y2)
+
+	    # Mark position on Y axis.
+	    if (abs(y) > .001)
+		call imp_markpos (gp, nint(y), ny)
+	} else {
+	    call glabax (gp, Memc[plot_title], xlabel, "")
+	    call ggwind (gp, x1, x2, y1, y2)
+	}
+
+	# Draw data vector.
+	call gpline (gp, xv, yv, nx)
+
+	call sfree (sp)
+end
+
+
+# IMP_SWIND -- Set all or part of the plotting window if autoscaling is not
+# desired.
+
+procedure imp_swind (n_x1, n_x2, n_y1, n_y2)
+
+real	n_x1, n_x2		# range of world coords in X
+real	n_y1, n_y2		# range of world coords in Y
+
+real	x1, x2, y1, y2
+common	/implcom/ x1, x2, y1, y2
+
+begin
+	if (!IS_INDEF(n_x1))
+	    x1 = n_x1
+	if (!IS_INDEF(n_x2))
+	    x2 = n_x2
+	if (!IS_INDEF(n_y1))
+	    y1 = n_y1
+	if (!IS_INDEF(n_y2))
+	    y2 = n_y2
+end
+
+
+# IMP_REDRAW -- Erase the old vector and draw a new one in its place.
+
+procedure imp_redraw (gp, xold, yold, xnew, ynew, npix)
+
+pointer	gp			# graphics descriptor
+real	xold[ARB], yold[ARB]	# old data vector
+real	xnew[ARB], ynew[ARB]	# new data vector
+int	npix			# length of the data vectors
+
+int	i, n
+
+begin
+	# Erase the old vector and redraw the new in its place, in segments
+	# of length SEGLEN.  These segments must overlap by one pixel to
+	# produce a continuous output polyline.
+
+	do i = 1, npix, SEGLEN {
+	    n  = min (SEGLEN + 1, npix - i + 1)
+
+	    # Erase next segment of old vector.
+	    call gseti (gp, G_PLTYPE, 0)
+	    call gpline (gp, xold[i], yold[i], n)
+
+	    # Plot same segment of new vector.
+	    call gseti (gp, G_PLTYPE, 1)
+	    call gseti (gp, G_PLCOLOR, 1)
+	    call gpline (gp, xnew[i], ynew[i], n)
+	}
+end
+
+
+# IMP_MARKPOS -- Mark the line or column number on the right axis of the plot.
+
+procedure imp_markpos (gp, line, nlines)
+
+pointer	gp			# graphics descriptor
+int	line			# line or column
+int	nlines			# number of lines or columns
+real	y, x1, x2, y1, y2
+
+begin
+	if (nlines < 2)
+	    return
+
+	call ggwind (gp, x1, x2, y1, y2)
+	y = (y2 - y1) / (nlines - 1) * (line - 1) + y1
+	call gmark (gp, x2, y, GM_PLUS, 3., 4.)
+end
+
+
+# IMP_TRAN -- Transform cursor coordinate to line and column in image.
+
+procedure imp_tran (gp, x, y, px, py, xvec, nx, ny)
+
+pointer	gp			# graphics descriptor
+real	x, y			# cursor coordinate
+real	px, py			# image coordinate
+real	xvec[nx]		# x vector
+int	nx, ny			# number of columns and lines
+
+int	i
+real	x1, x2, y1, y2, diff, diffmin
+bool	fp_equalr()
+
+begin
+	call ggwind (gp, x1, x2, y1, y2)
+	if (fp_equalr (y1, y2))
+	    py = nint (ny / 2.)
+	else
+	    py = nint ((ny - 1) / (y2 - y1) * (y - y1) + 1)
+
+	px = 1
+	diffmin = abs (x - xvec[1])
+	do i = 2, nx {
+	    diff = abs (x - xvec[i])
+	    if (diff < diffmin) {
+		px = i
+		diffmin = diff
+	    }
+	}
+end