Initial commit

37 files changed, 2483 insertions, 0 deletions

diff --git a/pkg/xtools/ranges/Revisions b/pkg/xtools/ranges/Revisions
new file mode 100644
index 00000000..a249d529
--- /dev/null
+++ b/pkg/xtools/ranges/Revisions
@@ -0,0 +1,59 @@
+.help revisions Jun88 pkg.xtools.ranges
+.nf
+xtools$ranges/rgranges.x
+    The range parsing would fail for single numbers.  (4/29/94, Valdes)
+
+xtools$ranges/rgranges.x
+xtools$ranges/rgxranges.gx
+    Added @file capability.  Rewrote the parsing logic.  (10/9/91, Valdes)
+
+xtools$ranges/rgwtbin.gx
+    At least one point from each sample region must be included regardless
+    of the size of average.  (6/14/89, Valdes)
+
+xtools$ranges/rgwtbin.gx
+    The remainder bin in a sample region must be at least
+    max (min (N, 3), (N+1)/2) except that a single bin may be any size. 
+    (6/7/89, Valdes)
+ 
+xtools$ranges/rgencode.x, rginverse.x, rgnext.x, rgunion.x
+    Added some missing functionality to convert a range into a string,
+    to invert a range, to get the next higher member of a range (ala
+    xtools$ranges.x) and to take the union of two ranges.
+    (6/2/89, Seaman)
+
+xtools$ranges/rgintersect.x, rgmerge.x
+    Fixed bugs in handling overlapping ranges.  (6/2/89, Seaman)
+
+xtools$ranges/rgxranges.gx
+    Numbers in scientific notation are now recognized.  Based on report
+    from Ivo Busko.  (3/1/89)
+ 
+xtools$ranges/rgwtbin.gx
+    The remainder bin in a sample region must be at least max (3, (N+1)/2)
+    except that a single bin may be any size.  (1/23/89, Valdes)
+ 
+xtools$ranges/rgdump.x +
+xtools$ranges/rgmerge.x
+    Valdes, May 4, 1987
+    1.  Added a debugging procedure for dumping the ranges descriptor.
+    2.  Fixed a bug when merging overlapping ranges.
+
+xtools$ranges/rgbin.gx
+xtools$ranges/rgwtbin.gx
+    Valdes, August 11, 1986
+    1.  Since AMED$T no longer modifies the input array the temporary arrays
+	used to preserve the input array are no longer needed.
+
+xtools$ranges: Valdes, August 11, 1986
+    1.  Reorganized package to have separate modules for each datatype.
+	This allows loading only the required procedures.
+
+xtools$ranges/rgwtbin.gx: Valdes, August 8, 1986
+    1.  If all the weights were zero in a given range then a divide by zero
+	would result.  A check against this was added.
+
+xtools$ranges: Valdes, March 13, 1986
+    1.  The RANGES package has been converted to generic form.  It is compiled
+	into both single and double precision procedures.
+.endhelp
diff --git a/pkg/xtools/ranges/mkpkg b/pkg/xtools/ranges/mkpkg
new file mode 100644
index 00000000..9ec1673f
--- /dev/null
+++ b/pkg/xtools/ranges/mkpkg
@@ -0,0 +1,49 @@
+# Library for the RANGES procedures.
+
+$checkout libxtools.a lib$
+$update	  libxtools.a
+$checkin  libxtools.a lib$
+$exit
+
+generic:
+	$set	GEN = "$$generic -k -t rd"
+	$ifolder (rgbinr.x, rgbin.gx)		$(GEN) rgbin.gx $endif
+	$ifolder (rgexcluder.x, rgexclude.gx)	$(GEN) rgexclude.gx $endif
+	$ifolder (rggxmarkr.x, rggxmark.gx)	$(GEN) rggxmark.gx $endif
+	$ifolder (rgpackr.x, rgpack.gx)		$(GEN) rgpack.gx $endif
+	$ifolder (rgunpackr.x, rgunpack.gx)	$(GEN) rgunpack.gx $endif
+	$ifolder (rgwtbinr.x, rgwtbin.gx)	$(GEN) rgwtbin.gx $endif
+	$ifolder (rgxrangesr.x, rgxranges.gx)	$(GEN) rgxranges.gx $endif
+	;
+
+libxtools.a:
+	$ifeq (USE_GENERIC, yes) $call generic $endif
+
+	rgbind.x	<pkg/rg.h>
+	rgbinr.x	<pkg/rg.h>
+	rgdump.x	<pkg/rg.h>
+	rgencode.x	<pkg/rg.h>
+	rgexcluded.x	<pkg/rg.h>
+	rgexcluder.x	<pkg/rg.h>
+	rgfree.x	
+	rggxmarkd.x	<gset.h> <pkg/rg.h>
+	rggxmarkr.x	<gset.h> <pkg/rg.h>
+	rgindices.x	<pkg/rg.h>
+	rginrange.x	<pkg/rg.h>
+	rgintersect.x	<pkg/rg.h>
+	rginverse.x	<pkg/rg.h>
+	rgmerge.x	<pkg/rg.h>
+	rgnext.x	<mach.h> <pkg/rg.h>
+	rgorder.x	<pkg/rg.h>
+	rgpackd.x	<pkg/rg.h>
+	rgpackr.x	<pkg/rg.h>
+	rgranges.x	<ctype.h> <error.h> <pkg/rg.h>
+	rgunion.x	<pkg/rg.h>
+	rgunpackd.x	<pkg/rg.h>
+	rgunpackr.x	<pkg/rg.h>
+	rgwindow.x	<pkg/rg.h>
+	rgwtbind.x	<pkg/rg.h>
+	rgwtbinr.x	<pkg/rg.h>
+	rgxrangesd.x	<ctype.h> <error.h> <pkg/rg.h>
+	rgxrangesr.x	<ctype.h> <error.h> <pkg/rg.h>
+	;
diff --git a/pkg/xtools/ranges/rgbin.gx b/pkg/xtools/ranges/rgbin.gx
new file mode 100644
index 00000000..f1133a1c
--- /dev/null
+++ b/pkg/xtools/ranges/rgbin.gx
@@ -0,0 +1,75 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<pkg/rg.h>
+
+# RG_BIN -- Average or median of data.
+#
+# The ranges are broken up into subranges of at most abs (nbin) points.  The
+# subranges are averaged if nbin > 1 and medianed if nbin < 1.
+# The output array must be large enough to contain the desired points.
+# If the ranges are merged then the input and output arrays may be the same.
+
+procedure rg_bin$t (rg, nbin, in, nin, out, nout)
+
+pointer	rg				# Ranges
+int	nbin				# Maximum points in average or median
+PIXEL	in[nin]				# Input array
+int	nin				# Number of input points
+PIXEL	out[ARB]			# Output array
+int	nout				# Number of output points
+
+int	i, j, k, n, npts, ntemp
+
+PIXEL	asum$t(), amed$t()
+
+errchk	rg_pack$t
+
+begin
+	# Error check the range pointer.
+
+	if (rg == NULL)
+	    call error (0, "Range descriptor undefined")
+
+	# If the bin size is exactly one then move the selected input points
+	# to the output array.
+
+	if (abs (nbin) == 1) {
+	    call rg_pack$t (rg, in, out)
+	    return
+	}
+
+	# Determine the subranges and take the median or average.
+
+	npts = abs (nbin)
+	ntemp = 0
+
+	do i = 1, RG_NRGS(rg) {
+	    if (RG_X1(rg, i) > RG_X2(rg, i)) {
+		j = min (nin, RG_X1(rg, i))
+		k = max (1, RG_X2(rg, i))
+		while (j >= k) {
+		    n = max (0, min (npts, j - k + 1))
+		    k = k - n
+		    ntemp = ntemp + 1
+		    if (nbin > 0)
+		        out[ntemp] = asum$t (in[k + 1], n) / n
+		    else
+			out[ntemp] = amed$t (in[k+1], n)
+		}
+	    } else {
+		j = max (1, RG_X1(rg, i))
+		k = min (nin, RG_X2(rg, i))
+	        while (j <= k) {
+		    n = max (0, min (npts, k - j + 1))
+		    ntemp = ntemp + 1
+		    if (nbin > 0)
+		        out[ntemp] = asum$t (in[j], n) / n
+		    else
+			out[ntemp] = amed$t (in[j], n)
+		    j = j + n
+		}
+	    }
+	}
+
+	nout = ntemp
+end
diff --git a/pkg/xtools/ranges/rgbind.x b/pkg/xtools/ranges/rgbind.x
new file mode 100644
index 00000000..16c66760
--- /dev/null
+++ b/pkg/xtools/ranges/rgbind.x
@@ -0,0 +1,75 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<pkg/rg.h>
+
+# RG_BIN -- Average or median of data.
+#
+# The ranges are broken up into subranges of at most abs (nbin) points.  The
+# subranges are averaged if nbin > 1 and medianed if nbin < 1.
+# The output array must be large enough to contain the desired points.
+# If the ranges are merged then the input and output arrays may be the same.
+
+procedure rg_bind (rg, nbin, in, nin, out, nout)
+
+pointer	rg				# Ranges
+int	nbin				# Maximum points in average or median
+double	in[nin]				# Input array
+int	nin				# Number of input points
+double	out[ARB]			# Output array
+int	nout				# Number of output points
+
+int	i, j, k, n, npts, ntemp
+
+double	asumd(), amedd()
+
+errchk	rg_packd
+
+begin
+	# Error check the range pointer.
+
+	if (rg == NULL)
+	    call error (0, "Range descriptor undefined")
+
+	# If the bin size is exactly one then move the selected input points
+	# to the output array.
+
+	if (abs (nbin) == 1) {
+	    call rg_packd (rg, in, out)
+	    return
+	}
+
+	# Determine the subranges and take the median or average.
+
+	npts = abs (nbin)
+	ntemp = 0
+
+	do i = 1, RG_NRGS(rg) {
+	    if (RG_X1(rg, i) > RG_X2(rg, i)) {
+		j = min (nin, RG_X1(rg, i))
+		k = max (1, RG_X2(rg, i))
+		while (j >= k) {
+		    n = max (0, min (npts, j - k + 1))
+		    k = k - n
+		    ntemp = ntemp + 1
+		    if (nbin > 0)
+		        out[ntemp] = asumd (in[k + 1], n) / n
+		    else
+			out[ntemp] = amedd (in[k+1], n)
+		}
+	    } else {
+		j = max (1, RG_X1(rg, i))
+		k = min (nin, RG_X2(rg, i))
+	        while (j <= k) {
+		    n = max (0, min (npts, k - j + 1))
+		    ntemp = ntemp + 1
+		    if (nbin > 0)
+		        out[ntemp] = asumd (in[j], n) / n
+		    else
+			out[ntemp] = amedd (in[j], n)
+		    j = j + n
+		}
+	    }
+	}
+
+	nout = ntemp
+end
diff --git a/pkg/xtools/ranges/rgbinr.x b/pkg/xtools/ranges/rgbinr.x
new file mode 100644
index 00000000..81fb9f70
--- /dev/null
+++ b/pkg/xtools/ranges/rgbinr.x
@@ -0,0 +1,75 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<pkg/rg.h>
+
+# RG_BIN -- Average or median of data.
+#
+# The ranges are broken up into subranges of at most abs (nbin) points.  The
+# subranges are averaged if nbin > 1 and medianed if nbin < 1.
+# The output array must be large enough to contain the desired points.
+# If the ranges are merged then the input and output arrays may be the same.
+
+procedure rg_binr (rg, nbin, in, nin, out, nout)
+
+pointer	rg				# Ranges
+int	nbin				# Maximum points in average or median
+real	in[nin]				# Input array
+int	nin				# Number of input points
+real	out[ARB]			# Output array
+int	nout				# Number of output points
+
+int	i, j, k, n, npts, ntemp
+
+real	asumr(), amedr()
+
+errchk	rg_packr
+
+begin
+	# Error check the range pointer.
+
+	if (rg == NULL)
+	    call error (0, "Range descriptor undefined")
+
+	# If the bin size is exactly one then move the selected input points
+	# to the output array.
+
+	if (abs (nbin) == 1) {
+	    call rg_packr (rg, in, out)
+	    return
+	}
+
+	# Determine the subranges and take the median or average.
+
+	npts = abs (nbin)
+	ntemp = 0
+
+	do i = 1, RG_NRGS(rg) {
+	    if (RG_X1(rg, i) > RG_X2(rg, i)) {
+		j = min (nin, RG_X1(rg, i))
+		k = max (1, RG_X2(rg, i))
+		while (j >= k) {
+		    n = max (0, min (npts, j - k + 1))
+		    k = k - n
+		    ntemp = ntemp + 1
+		    if (nbin > 0)
+		        out[ntemp] = asumr (in[k + 1], n) / n
+		    else
+			out[ntemp] = amedr (in[k+1], n)
+		}
+	    } else {
+		j = max (1, RG_X1(rg, i))
+		k = min (nin, RG_X2(rg, i))
+	        while (j <= k) {
+		    n = max (0, min (npts, k - j + 1))
+		    ntemp = ntemp + 1
+		    if (nbin > 0)
+		        out[ntemp] = asumr (in[j], n) / n
+		    else
+			out[ntemp] = amedr (in[j], n)
+		    j = j + n
+		}
+	    }
+	}
+
+	nout = ntemp
+end
diff --git a/pkg/xtools/ranges/rgdump.x b/pkg/xtools/ranges/rgdump.x
new file mode 100644
index 00000000..97c3a89b
--- /dev/null
+++ b/pkg/xtools/ranges/rgdump.x
@@ -0,0 +1,28 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+
+include	<pkg/rg.h>
+
+# RG_DUMP -- Dump the contents of a range structure.
+
+procedure rg_dump (rg)
+
+pointer	rg		# Ranges
+
+int	i
+
+begin
+	if (rg == NULL)
+	    call printf ("RG_DUMP: The range pointer is NULL\n")
+	else {
+	    call printf ("RG_DUMP: NPTS = %d, NRGS = %d\n")
+		call pargi (RG_NPTS(rg))
+		call pargi (RG_NRGS(rg))
+	    do i = 1, RG_NRGS(rg) {
+		call printf ("  %4d - %4d\n")
+		    call pargi (RG_X1(rg, i))
+		    call pargi (RG_X2(rg, i))
+	    }
+	}
+	call flush (STDOUT)
+end
diff --git a/pkg/xtools/ranges/rgencode.x b/pkg/xtools/ranges/rgencode.x
new file mode 100644
index 00000000..ff0a0343
--- /dev/null
+++ b/pkg/xtools/ranges/rgencode.x
@@ -0,0 +1,52 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<pkg/rg.h>
+
+# RG_ENCODE -- Encode a range structure into a string, return the
+# number of characters that were written or ERR for string overflow.
+
+int procedure rg_encode (rg, outstr, maxch)
+
+pointer	rg				# First set of ranges
+char	outstr[maxch]			# String to receive the ranges
+int	maxch				# Maximum length of the string
+
+char	tmpstr[SZ_LINE]
+int	i, outlen
+
+int	strlen()
+
+begin
+	if (rg == NULL)
+	    call error (0, "Range descriptor undefined")
+
+	outlen = 0
+	outstr[1] = EOS
+
+	do i = 1, RG_NRGS(rg) {
+	    if (RG_X1(rg, i) != RG_X2(rg, i)) {
+		call sprintf (tmpstr, maxch, "%d:%d,")
+		    call pargi (RG_X1(rg, i))
+		    call pargi (RG_X2(rg, i))
+	    } else {
+		call sprintf (tmpstr, maxch, "%d,")
+		    call pargi (RG_X1(rg, i))
+	    }
+
+	    outlen = outlen + strlen (tmpstr)
+
+	    if (outlen <= maxch)
+		call strcat (tmpstr, outstr, maxch)
+	    else {
+		outstr[1] = EOS
+		return (ERR)
+	    }
+	}
+
+	# remove the last comma
+
+	outstr[outlen] = EOS
+	outlen = outlen - 1
+
+	return (outlen)
+end
diff --git a/pkg/xtools/ranges/rgexclude.gx b/pkg/xtools/ranges/rgexclude.gx
new file mode 100644
index 00000000..876e4ef7
--- /dev/null
+++ b/pkg/xtools/ranges/rgexclude.gx
@@ -0,0 +1,56 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<pkg/rg.h>
+
+# RG_EXCLUDE -- Exclude points given by ranges.
+#
+# The output array must be large enough to contain the desired points.
+# If the ranges are merged then the input and output arrays may be the same.
+
+procedure rg_exclude$t (rg, a, nin, b, nout)
+
+pointer	rg					# Ranges
+PIXEL	a[nin]					# Input array
+int	nin					# Number of input points
+PIXEL	b[ARB]					# Output array
+int	nout					# Number of output points
+
+int	i, j, k, n, ntemp
+
+begin
+	# Error check the range pointer.
+
+	if (rg == NULL)
+	    call error (0, "Range descriptor undefined")
+
+	if (RG_NRGS(rg) == 0) {
+	    call amov$t (a[1], b[1], nin)
+	    nout = nin
+	} else {
+	    ntemp = 0
+
+	    i = 1
+	    j = 1
+	    k = min (nin, min (RG_X1(rg, i), RG_X2(rg, i)) - 1)
+	    n = max (0, k - j + 1)
+	    call amov$t (a[j], b[ntemp+1], n)
+	    ntemp = ntemp + n
+
+	    do i = 2, RG_NRGS(rg) {
+	        j = max (1, max (RG_X1(rg, i-1), RG_X2(rg, i-1)) + 1)
+	        k = min (nin, min (RG_X1(rg, i), RG_X2(rg, i)) - 1)
+	        n = max (0, k - j + 1)
+	        call amov$t (a[j], b[ntemp+1], n)
+	        ntemp = ntemp + n
+	    }
+
+	    i = RG_NRGS (rg)
+	    j = max (1, max (RG_X1(rg, i), RG_X2(rg, i)) + 1)
+	    k = nin
+	    n = max (0, k - j + 1)
+	    call amov$t (a[j], b[ntemp+1], n)
+	    ntemp = ntemp + n
+	}
+
+	nout = ntemp
+end
diff --git a/pkg/xtools/ranges/rgexcluded.x b/pkg/xtools/ranges/rgexcluded.x
new file mode 100644
index 00000000..2d9ef823
--- /dev/null
+++ b/pkg/xtools/ranges/rgexcluded.x
@@ -0,0 +1,56 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<pkg/rg.h>
+
+# RG_EXCLUDE -- Exclude points given by ranges.
+#
+# The output array must be large enough to contain the desired points.
+# If the ranges are merged then the input and output arrays may be the same.
+
+procedure rg_excluded (rg, a, nin, b, nout)
+
+pointer	rg					# Ranges
+double	a[nin]					# Input array
+int	nin					# Number of input points
+double	b[ARB]					# Output array
+int	nout					# Number of output points
+
+int	i, j, k, n, ntemp
+
+begin
+	# Error check the range pointer.
+
+	if (rg == NULL)
+	    call error (0, "Range descriptor undefined")
+
+	if (RG_NRGS(rg) == 0) {
+	    call amovd (a[1], b[1], nin)
+	    nout = nin
+	} else {
+	    ntemp = 0
+
+	    i = 1
+	    j = 1
+	    k = min (nin, min (RG_X1(rg, i), RG_X2(rg, i)) - 1)
+	    n = max (0, k - j + 1)
+	    call amovd (a[j], b[ntemp+1], n)
+	    ntemp = ntemp + n
+
+	    do i = 2, RG_NRGS(rg) {
+	        j = max (1, max (RG_X1(rg, i-1), RG_X2(rg, i-1)) + 1)
+	        k = min (nin, min (RG_X1(rg, i), RG_X2(rg, i)) - 1)
+	        n = max (0, k - j + 1)
+	        call amovd (a[j], b[ntemp+1], n)
+	        ntemp = ntemp + n
+	    }
+
+	    i = RG_NRGS (rg)
+	    j = max (1, max (RG_X1(rg, i), RG_X2(rg, i)) + 1)
+	    k = nin
+	    n = max (0, k - j + 1)
+	    call amovd (a[j], b[ntemp+1], n)
+	    ntemp = ntemp + n
+	}
+
+	nout = ntemp
+end
diff --git a/pkg/xtools/ranges/rgexcluder.x b/pkg/xtools/ranges/rgexcluder.x
new file mode 100644
index 00000000..44cb90fe
--- /dev/null
+++ b/pkg/xtools/ranges/rgexcluder.x
@@ -0,0 +1,56 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<pkg/rg.h>
+
+# RG_EXCLUDE -- Exclude points given by ranges.
+#
+# The output array must be large enough to contain the desired points.
+# If the ranges are merged then the input and output arrays may be the same.
+
+procedure rg_excluder (rg, a, nin, b, nout)
+
+pointer	rg					# Ranges
+real	a[nin]					# Input array
+int	nin					# Number of input points
+real	b[ARB]					# Output array
+int	nout					# Number of output points
+
+int	i, j, k, n, ntemp
+
+begin
+	# Error check the range pointer.
+
+	if (rg == NULL)
+	    call error (0, "Range descriptor undefined")
+
+	if (RG_NRGS(rg) == 0) {
+	    call amovr (a[1], b[1], nin)
+	    nout = nin
+	} else {
+	    ntemp = 0
+
+	    i = 1
+	    j = 1
+	    k = min (nin, min (RG_X1(rg, i), RG_X2(rg, i)) - 1)
+	    n = max (0, k - j + 1)
+	    call amovr (a[j], b[ntemp+1], n)
+	    ntemp = ntemp + n
+
+	    do i = 2, RG_NRGS(rg) {
+	        j = max (1, max (RG_X1(rg, i-1), RG_X2(rg, i-1)) + 1)
+	        k = min (nin, min (RG_X1(rg, i), RG_X2(rg, i)) - 1)
+	        n = max (0, k - j + 1)
+	        call amovr (a[j], b[ntemp+1], n)
+	        ntemp = ntemp + n
+	    }
+
+	    i = RG_NRGS (rg)
+	    j = max (1, max (RG_X1(rg, i), RG_X2(rg, i)) + 1)
+	    k = nin
+	    n = max (0, k - j + 1)
+	    call amovr (a[j], b[ntemp+1], n)
+	    ntemp = ntemp + n
+	}
+
+	nout = ntemp
+end
diff --git a/pkg/xtools/ranges/rgfree.x b/pkg/xtools/ranges/rgfree.x
new file mode 100644
index 00000000..8b2ab344
--- /dev/null
+++ b/pkg/xtools/ranges/rgfree.x
@@ -0,0 +1,14 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+# RG_FREE -- Free ranges memory.
+
+procedure rg_free (rg)
+
+pointer	rg				# Ranges
+
+begin
+	if (rg != NULL) {
+	    call mfree (rg, TY_STRUCT)
+	    rg = NULL
+	}
+end
diff --git a/pkg/xtools/ranges/rggxmark.gx b/pkg/xtools/ranges/rggxmark.gx
new file mode 100644
index 00000000..26108c98
--- /dev/null
+++ b/pkg/xtools/ranges/rggxmark.gx
@@ -0,0 +1,52 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<gset.h>
+include	<pkg/rg.h>
+
+# RG_GXMARK -- Mark x ranges.
+
+procedure rg_gxmark$t (gp, rstr, x, npts, pltype)
+
+pointer	gp				# GIO pointer
+char	rstr[ARB]			# Range string
+PIXEL	x[npts]				# Ordinates of graph
+int	npts				# Number of data points
+int	pltype				# Plot line type
+
+pointer	rg
+int	i, pltype1
+real	xl, xr, yb, yt, dy
+real	x1, x2, y1, y2, y3
+
+int	gstati(), stridxs()
+pointer	rg_xranges$t()
+
+begin
+	if (stridxs ("*", rstr) > 0)
+	    return
+
+	rg = rg_xranges$t (rstr, x, npts)
+
+	pltype1 = gstati (gp, G_PLTYPE)
+	call gseti (gp, G_PLTYPE, pltype)
+
+	call ggwind (gp, xl, xr, yb, yt)
+
+	dy = yt - yb
+	y1 = yb + dy / 100
+	y2 = y1 + dy / 20
+	y3 = (y1 + y2) / 2
+
+	do i = 1, RG_NRGS(rg) {
+	    x1 = x[RG_X1(rg, i)]
+	    x2 = x[RG_X2(rg, i)]
+	    if ((x1 > xl) && (x1 < xr))
+	        call gline (gp, x1, y1, x1, y2)
+	    if ((x2 > xl) && (x2 < xr))
+	        call gline (gp, x2, y1, x2, y2)
+	    call gline (gp, x1, y3, x2, y3)
+	}
+
+	call gseti (gp, G_PLTYPE, pltype1)
+	call rg_free (rg)
+end
diff --git a/pkg/xtools/ranges/rggxmarkd.x b/pkg/xtools/ranges/rggxmarkd.x
new file mode 100644
index 00000000..82eb49db
--- /dev/null
+++ b/pkg/xtools/ranges/rggxmarkd.x
@@ -0,0 +1,52 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<gset.h>
+include	<pkg/rg.h>
+
+# RG_GXMARK -- Mark x ranges.
+
+procedure rg_gxmarkd (gp, rstr, x, npts, pltype)
+
+pointer	gp				# GIO pointer
+char	rstr[ARB]			# Range string
+double	x[npts]				# Ordinates of graph
+int	npts				# Number of data points
+int	pltype				# Plot line type
+
+pointer	rg
+int	i, pltype1
+real	xl, xr, yb, yt, dy
+real	x1, x2, y1, y2, y3
+
+int	gstati(), stridxs()
+pointer	rg_xrangesd()
+
+begin
+	if (stridxs ("*", rstr) > 0)
+	    return
+
+	rg = rg_xrangesd (rstr, x, npts)
+
+	pltype1 = gstati (gp, G_PLTYPE)
+	call gseti (gp, G_PLTYPE, pltype)
+
+	call ggwind (gp, xl, xr, yb, yt)
+
+	dy = yt - yb
+	y1 = yb + dy / 100
+	y2 = y1 + dy / 20
+	y3 = (y1 + y2) / 2
+
+	do i = 1, RG_NRGS(rg) {
+	    x1 = x[RG_X1(rg, i)]
+	    x2 = x[RG_X2(rg, i)]
+	    if ((x1 > xl) && (x1 < xr))
+	        call gline (gp, x1, y1, x1, y2)
+	    if ((x2 > xl) && (x2 < xr))
+	        call gline (gp, x2, y1, x2, y2)
+	    call gline (gp, x1, y3, x2, y3)
+	}
+
+	call gseti (gp, G_PLTYPE, pltype1)
+	call rg_free (rg)
+end
diff --git a/pkg/xtools/ranges/rggxmarkr.x b/pkg/xtools/ranges/rggxmarkr.x
new file mode 100644
index 00000000..ec0f63b8
--- /dev/null
+++ b/pkg/xtools/ranges/rggxmarkr.x
@@ -0,0 +1,52 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<gset.h>
+include	<pkg/rg.h>
+
+# RG_GXMARK -- Mark x ranges.
+
+procedure rg_gxmarkr (gp, rstr, x, npts, pltype)
+
+pointer	gp				# GIO pointer
+char	rstr[ARB]			# Range string
+real	x[npts]				# Ordinates of graph
+int	npts				# Number of data points
+int	pltype				# Plot line type
+
+pointer	rg
+int	i, pltype1
+real	xl, xr, yb, yt, dy
+real	x1, x2, y1, y2, y3
+
+int	gstati(), stridxs()
+pointer	rg_xrangesr()
+
+begin
+	if (stridxs ("*", rstr) > 0)
+	    return
+
+	rg = rg_xrangesr (rstr, x, npts)
+
+	pltype1 = gstati (gp, G_PLTYPE)
+	call gseti (gp, G_PLTYPE, pltype)
+
+	call ggwind (gp, xl, xr, yb, yt)
+
+	dy = yt - yb
+	y1 = yb + dy / 100
+	y2 = y1 + dy / 20
+	y3 = (y1 + y2) / 2
+
+	do i = 1, RG_NRGS(rg) {
+	    x1 = x[RG_X1(rg, i)]
+	    x2 = x[RG_X2(rg, i)]
+	    if ((x1 > xl) && (x1 < xr))
+	        call gline (gp, x1, y1, x1, y2)
+	    if ((x2 > xl) && (x2 < xr))
+	        call gline (gp, x2, y1, x2, y2)
+	    call gline (gp, x1, y3, x2, y3)
+	}
+
+	call gseti (gp, G_PLTYPE, pltype1)
+	call rg_free (rg)
+end
diff --git a/pkg/xtools/ranges/rgindices.x b/pkg/xtools/ranges/rgindices.x
new file mode 100644
index 00000000..48f1ec8f
--- /dev/null
+++ b/pkg/xtools/ranges/rgindices.x
@@ -0,0 +1,81 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<pkg/rg.h>
+
+# RG_INDICES -- Return the indices in the ranges.
+
+procedure rg_indices (rg, indices, npts, type)
+
+pointer	rg				# Ranges
+pointer	indices				# Indices
+int	npts				# Number of indices
+int	type				# Data type of points
+
+int	i, j, k, step
+
+begin
+	# Error check the range pointer.
+
+	if (rg == NULL)
+	    call error (0, "Range descriptor undefined")
+
+	# Determine the number of range points.
+
+	indices = NULL
+	npts = 0
+	if (RG_NRGS (rg) == 0)
+	    return
+
+	do i = 1, RG_NRGS(rg) {
+	    if (RG_X1(rg, i) > RG_X2(rg, i))
+		npts = npts + RG_X1(rg, i) - RG_X2(rg, i) + 1
+	    else
+		npts = npts + RG_X2(rg, i) - RG_X1(rg, i) + 1
+	}
+
+	# Allocate the range points array.
+
+	call malloc (indices, npts, type)
+
+	# Set the range points.
+
+	k = indices
+	do i = 1, RG_NRGS(rg) {
+	    if (RG_X1(rg, i) > RG_X2(rg, i))
+		step = -1
+	    else
+		step = 1
+
+	    switch (type) {
+	    case TY_SHORT:
+	    	do j = RG_X1(rg, i), RG_X2(rg, i), step {
+		    Mems[k] = j
+		    k = k + 1
+		}
+	    case TY_INT:
+	    	do j = RG_X1(rg, i), RG_X2(rg, i), step {
+		    Memi[k] = j
+		    k = k + 1
+		}
+	    case TY_LONG:
+	    	do j = RG_X1(rg, i), RG_X2(rg, i), step {
+		    Meml[k] = j
+		    k = k + 1
+		}
+	    case TY_REAL:
+	    	do j = RG_X1(rg, i), RG_X2(rg, i), step {
+		    Memr[k] = j
+		    k = k + 1
+		}
+	    case TY_DOUBLE:
+	    	do j = RG_X1(rg, i), RG_X2(rg, i), step {
+		    Memd[k] = j
+		    k = k + 1
+		}
+	    default:
+		call error (0, "rg_indices:  Datatype not available")
+	    }
+	}
+
+	return
+end
diff --git a/pkg/xtools/ranges/rginrange.x b/pkg/xtools/ranges/rginrange.x
new file mode 100644
index 00000000..7dd946ae
--- /dev/null
+++ b/pkg/xtools/ranges/rginrange.x
@@ -0,0 +1,29 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<pkg/rg.h>
+
+# RG_INRANGE -- Is value in the ranges?
+
+int procedure rg_inrange (rg, rval)
+
+pointer	rg				# Ranges
+int	rval				# Range value to test
+
+int	i
+
+begin
+	# Error check the range pointer.
+
+	if (rg == NULL)
+	    call error (0, "Range descriptor undefined")
+
+	do i = 1, RG_NRGS(rg) {
+	    if ((RG_X1(rg, i) <= RG_X2(rg, i)) && (rval >= RG_X1(rg, i)) &&
+		    (rval <= RG_X2(rg, i)))
+		return (YES)
+	    else if ((rval >= RG_X2(rg, i)) && (rval <= RG_X1(rg, i)))
+		return (YES)
+	}
+
+	return (NO)
+end
diff --git a/pkg/xtools/ranges/rgintersect.x b/pkg/xtools/ranges/rgintersect.x
new file mode 100644
index 00000000..5e4e4390
--- /dev/null
+++ b/pkg/xtools/ranges/rgintersect.x
@@ -0,0 +1,58 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<pkg/rg.h>
+
+# RG_INTERSECT -- Intersect two sets of ordered and merged ranges.
+
+pointer procedure rg_intersect (rg1, rg2)
+
+pointer	rg1				# First set of ranges
+pointer	rg2				# Second set of ranges
+
+pointer	rg3				# Pointer to intersection
+
+int	i, j, k
+
+begin
+	# Error check the range pointers.
+
+	if ((rg1 == NULL) || (rg2 == NULL))
+	    call error (0, "Range descriptor(s) undefined")
+
+	# Allocate the range points array.
+
+	k = RG_NRGS(rg1) + RG_NRGS(rg2) - 1
+	call malloc (rg3, LEN_RG + 2 * max (1, k), TY_STRUCT)
+
+	# Set the ranges.
+
+	i = 1
+	j = 1
+	k = 0
+
+	while (i <= RG_NRGS(rg1) && j <= RG_NRGS(rg2)) {
+	    if (RG_X2(rg1, i) < RG_X1(rg2, j))
+		i = i + 1
+	    else if (RG_X2(rg2, j) < RG_X1(rg1, i))
+		j = j + 1
+	    else {
+		k = k + 1
+		RG_X1(rg3, k) = max (RG_X1(rg1, i), RG_X1(rg2, j))
+		RG_X2(rg3, k) = min (RG_X2(rg1, i), RG_X2(rg2, j))
+
+	        if (RG_X2(rg1, i) < RG_X2(rg2, j))
+		    i = i + 1
+	        else
+		    j = j + 1
+	    }
+	}
+
+	call realloc (rg3, LEN_RG + 2 * max (1, k), TY_STRUCT)
+
+	RG_NRGS(rg3) = k
+	RG_NPTS(rg3) = 0
+	do i = 1, RG_NRGS(rg3)
+	    RG_NPTS(rg3) = RG_NPTS(rg3) + RG_X2(rg3, i) - RG_X1(rg3, i) + 1
+
+	return (rg3)
+end
diff --git a/pkg/xtools/ranges/rginverse.x b/pkg/xtools/ranges/rginverse.x
new file mode 100644
index 00000000..869fde19
--- /dev/null
+++ b/pkg/xtools/ranges/rginverse.x
@@ -0,0 +1,34 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<pkg/rg.h>
+
+# RG_INVERSE -- Invert a set of ordered and merged ranges.
+
+procedure rg_inverse (rg, rmin, rmax)
+
+pointer	rg			# RANGES pointer
+int	rmin			# Minimum value of window
+int	rmax			# Maximum value of window
+
+int	i
+pointer rgtmp
+
+pointer	rg_window()
+
+begin
+	call malloc (rgtmp, LEN_RG + 2 * (RG_NRGS(rg) + 1), TY_STRUCT)
+	RG_NRGS(rgtmp) = RG_NRGS(rg) + 1
+
+	RG_X1(rgtmp, 1) = rmin
+
+	do i = 1, RG_NRGS(rg) {
+	    RG_X2(rgtmp, i) = RG_X1(rg, i) - 1
+	    RG_X1(rgtmp, i+1) = RG_X2(rg, i) + 1
+	}
+
+	RG_X2(rgtmp, RG_NRGS(rgtmp)) = rmax
+
+	call rg_free (rg)
+	rg = rg_window (rgtmp, rmin, rmax)
+	call rg_free (rgtmp)
+end
diff --git a/pkg/xtools/ranges/rgmerge.x b/pkg/xtools/ranges/rgmerge.x
new file mode 100644
index 00000000..2cb5034a
--- /dev/null
+++ b/pkg/xtools/ranges/rgmerge.x
@@ -0,0 +1,38 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<pkg/rg.h>
+
+# RG_MERGE -- Merge overlapping ranges in set of ordered ranges.
+
+procedure rg_merge (rg)
+
+pointer	rg				# Ranges
+
+int	new, old
+
+begin
+	# Error check the range pointer.
+
+	if (rg == NULL)
+	    call error (0, "Range descriptor undefined")
+	if (RG_NRGS(rg) == 0)
+	    return
+
+	# Eliminate overlapping ranges and count the number of new ranges.
+
+	new = 1
+	do old = 2, RG_NRGS(rg)
+	    if (RG_X1(rg, old) > RG_X2(rg, new) + 1) {
+		new = new + 1
+		RG_X1(rg, new) = RG_X1(rg, old)
+		RG_X2(rg, new) = RG_X2(rg, old)
+	    } else
+		RG_X2(rg, new) = max (RG_X2(rg, old), RG_X2(rg, new))
+
+	call realloc (rg, LEN_RG + 2 * new, TY_STRUCT)
+
+	RG_NPTS(rg) = 0
+	RG_NRGS(rg) = new
+	do new = 1, RG_NRGS(rg)
+	    RG_NPTS(rg) = RG_NPTS(rg) + RG_X2(rg, new) - RG_X1(rg, new) + 1
+end
diff --git a/pkg/xtools/ranges/rgnext.x b/pkg/xtools/ranges/rgnext.x
new file mode 100644
index 00000000..354ef813
--- /dev/null
+++ b/pkg/xtools/ranges/rgnext.x
@@ -0,0 +1,32 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<mach.h>
+include	<pkg/rg.h>
+
+# RG_NEXT -- Return the next value in a set of ordered and merged ranges.
+# Number is set to the next value in the set of ranges or is unchanged
+# (and EOF is returned) if there are no more values.
+
+int procedure rg_next (rg, number)
+
+pointer	rg			# RANGES pointer
+int	number			# Both input and output parameter
+
+int	next_number, i
+
+begin
+	next_number = number + 1
+
+	do i = 1, RG_NRGS(rg)
+	    if (next_number > RG_X2(rg, i)) {
+		next
+	    } else if (next_number < RG_X1(rg, i)) {
+		number = RG_X1(rg, i)
+	    	return (number)
+	    } else {
+		number = next_number
+		return (number)
+	    }
+
+	return (EOF)
+end
diff --git a/pkg/xtools/ranges/rgorder.x b/pkg/xtools/ranges/rgorder.x
new file mode 100644
index 00000000..7864ecb2
--- /dev/null
+++ b/pkg/xtools/ranges/rgorder.x
@@ -0,0 +1,43 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include <pkg/rg.h>
+
+# RG_ORDER -- Make all ranges increasing and order by the starting point.
+
+procedure rg_order (rg)
+
+pointer	rg				# Ranges
+
+int	i, j, temp
+
+begin
+	# Error check the range pointer.
+
+	if (rg == NULL)
+	    call error (0, "Range descriptor undefined")
+
+	# Make all ranges increasing.
+
+	do i = 1, RG_NRGS(rg) {
+	    if (RG_X1(rg, i) > RG_X2(rg, i)) {
+		temp = RG_X1(rg, i)
+		RG_X1(rg, i) = RG_X2(rg, i)
+		RG_X2(rg, i) = temp
+	    }
+	}
+
+	# Sort the ranges in increasing order.
+
+	do i = 1, RG_NRGS(rg) - 1 {
+	    do j = i + 1, RG_NRGS(rg) {
+		if (RG_X1(rg, i) > RG_X1(rg, j)) {
+		    temp = RG_X1(rg, i)
+		    RG_X1(rg, i) = RG_X1(rg, j)
+		    RG_X1(rg, j) = temp
+		    temp = RG_X2(rg, i)
+		    RG_X2(rg, i) = RG_X2(rg, j)
+		    RG_X2(rg, j) = temp
+		}
+	    }
+	}
+end
diff --git a/pkg/xtools/ranges/rgpack.gx b/pkg/xtools/ranges/rgpack.gx
new file mode 100644
index 00000000..d77a3a09
--- /dev/null
+++ b/pkg/xtools/ranges/rgpack.gx
@@ -0,0 +1,37 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<pkg/rg.h>
+
+# RG_PACK -- Pack input data to include only points in the ranges.
+#
+# The output array must be large enough to contain the desired points.
+# If the ranges are merged then the input and output arrays may be the same.
+
+procedure rg_pack$t (rg, a, b)
+
+pointer	rg					# Ranges
+PIXEL	a[ARB]					# Input array
+PIXEL	b[ARB]					# Output array
+
+int	i, j, k, n
+
+begin
+	# Error check the range pointer.
+
+	if (rg == NULL)
+	    call error (0, "Range pointer undefined")
+
+	j = 0
+	do i = 1, RG_NRGS(rg) {
+	    if (RG_X1(rg, i) > RG_X2(rg, i)) {
+		do k = RG_X1(rg, i), RG_X2(rg, i), -1 {
+		    j = j + 1
+		    b[j] = a[k]
+		}
+	    } else {
+		n = RG_X2(rg, i) - RG_X1(rg, i) + 1
+	        call amov$t (a[RG_X1(rg, i)], b[j + 1], n)
+	        j = j + n
+	    }
+	}
+end
diff --git a/pkg/xtools/ranges/rgpackd.x b/pkg/xtools/ranges/rgpackd.x
new file mode 100644
index 00000000..a0889ec6
--- /dev/null
+++ b/pkg/xtools/ranges/rgpackd.x
@@ -0,0 +1,37 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<pkg/rg.h>
+
+# RG_PACK -- Pack input data to include only points in the ranges.
+#
+# The output array must be large enough to contain the desired points.
+# If the ranges are merged then the input and output arrays may be the same.
+
+procedure rg_packd (rg, a, b)
+
+pointer	rg					# Ranges
+double	a[ARB]					# Input array
+double	b[ARB]					# Output array
+
+int	i, j, k, n
+
+begin
+	# Error check the range pointer.
+
+	if (rg == NULL)
+	    call error (0, "Range pointer undefined")
+
+	j = 0
+	do i = 1, RG_NRGS(rg) {
+	    if (RG_X1(rg, i) > RG_X2(rg, i)) {
+		do k = RG_X1(rg, i), RG_X2(rg, i), -1 {
+		    j = j + 1
+		    b[j] = a[k]
+		}
+	    } else {
+		n = RG_X2(rg, i) - RG_X1(rg, i) + 1
+	        call amovd (a[RG_X1(rg, i)], b[j + 1], n)
+	        j = j + n
+	    }
+	}
+end
diff --git a/pkg/xtools/ranges/rgpackr.x b/pkg/xtools/ranges/rgpackr.x
new file mode 100644
index 00000000..e01307a3
--- /dev/null
+++ b/pkg/xtools/ranges/rgpackr.x
@@ -0,0 +1,37 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<pkg/rg.h>
+
+# RG_PACK -- Pack input data to include only points in the ranges.
+#
+# The output array must be large enough to contain the desired points.
+# If the ranges are merged then the input and output arrays may be the same.
+
+procedure rg_packr (rg, a, b)
+
+pointer	rg					# Ranges
+real	a[ARB]					# Input array
+real	b[ARB]					# Output array
+
+int	i, j, k, n
+
+begin
+	# Error check the range pointer.
+
+	if (rg == NULL)
+	    call error (0, "Range pointer undefined")
+
+	j = 0
+	do i = 1, RG_NRGS(rg) {
+	    if (RG_X1(rg, i) > RG_X2(rg, i)) {
+		do k = RG_X1(rg, i), RG_X2(rg, i), -1 {
+		    j = j + 1
+		    b[j] = a[k]
+		}
+	    } else {
+		n = RG_X2(rg, i) - RG_X1(rg, i) + 1
+	        call amovr (a[RG_X1(rg, i)], b[j + 1], n)
+	        j = j + n
+	    }
+	}
+end
diff --git a/pkg/xtools/ranges/rgranges.x b/pkg/xtools/ranges/rgranges.x
new file mode 100644
index 00000000..913fc2b9
--- /dev/null
+++ b/pkg/xtools/ranges/rgranges.x
@@ -0,0 +1,136 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<error.h>
+include	<ctype.h>
+include	<pkg/rg.h>
+
+define	NRGS	10		# Allocation size
+
+# RG_RANGES -- Parse a range string.  Return a pointer to the ranges.
+
+pointer procedure rg_ranges (rstr, rmin, rmax)
+
+char	rstr[ARB]		# Range string
+int	rmin			# Minimum value
+int	rmax			# Maximum value
+pointer	rg			# Range pointer
+
+int	i, fd, strlen(), open(), getline()
+pointer	sp, str, ptr
+errchk	open, rg_add
+
+begin
+	call smark (sp)
+	call salloc (str, max (strlen (rstr), SZ_LINE), TY_CHAR)
+	call calloc (rg, LEN_RG, TY_STRUCT)
+
+	i = 1
+	while (rstr[i] != EOS) {
+
+	    # Find beginning and end of a range and copy it to the work string
+	    while (IS_WHITE(rstr[i]) || rstr[i]==',' || rstr[i]=='\n')
+	        i = i + 1
+	    if (rstr[i] == EOS)
+		break
+
+	    ptr = str
+	    while (!(IS_WHITE(rstr[i]) || rstr[i]==',' || rstr[i]=='\n' ||
+		rstr[i]==EOS)) {
+		Memc[ptr] = rstr[i]
+	        i = i + 1
+		ptr = ptr + 1
+	    }
+	    Memc[ptr] = EOS
+
+	    # Add range(s)
+	    iferr {
+		if (Memc[str] == '@') {
+		    fd = open (Memc[str+1], READ_ONLY, TEXT_FILE)
+		    while (getline (fd, Memc[str]) != EOF) {
+			iferr (call rg_add (rg, Memc[str], rmin, rmax))
+			    call erract (EA_WARN)
+		    }
+		    call close (fd)
+		} else
+		    call rg_add (rg, Memc[str], rmin, rmax)
+	    } then
+		call erract (EA_WARN)
+	}
+
+	call sfree (sp)
+	return (rg)
+end
+
+
+# RG_ADD -- Add a range
+
+procedure rg_add (rg, rstr, rmin, rmax)
+
+pointer	rg			# Range descriptor
+char	rstr[ARB]		# Range string
+int	rmin			# Minimum value
+int	rmax			# Maximum value
+
+int	i, j, nrgs, strlen(), ctoi()
+int	rval1, rval2
+pointer	sp, str, ptr
+
+begin
+	call smark (sp)
+	call salloc (str, strlen (rstr), TY_CHAR)
+
+	i = 1
+	while (rstr[i] != EOS) {
+
+	    # Find beginning and end of a range and copy it to the work string
+	    while (IS_WHITE(rstr[i]) || rstr[i]==',' || rstr[i]=='\n')
+	        i = i + 1
+	    if (rstr[i] == EOS)
+		break
+
+	    ptr = str
+	    while (!(IS_WHITE(rstr[i]) || rstr[i]==',' || rstr[i]=='\n' ||
+		rstr[i]==EOS)) {
+		if (rstr[i] == ':')
+		    Memc[ptr] = ' '
+		else
+		    Memc[ptr] = rstr[i]
+	        i = i + 1
+		ptr = ptr + 1
+	    }
+	    Memc[ptr] = EOS
+
+	    # Parse range
+	    if (Memc[str] == '@')
+		call error (1, "Cannot nest @files")
+	    else if (Memc[str] == '*') {
+		rval1 = rmin
+		rval2 = rmax
+	    } else {
+		# Get range
+		j = 1
+		if (ctoi (Memc[str], j, rval1) == 0)
+		    call error (1, "Range syntax error")
+		if (ctoi (Memc[str], j, rval2) == 0)
+		    rval2 = rval1
+	    }
+
+	    # Check limits.
+	    j = rval1
+	    rval1 = min (j, rval2)
+	    rval2 = max (j, rval2)
+	    if (rval2 >= rmin && rval1 <= rmax) {
+		nrgs = RG_NRGS(rg)
+		if (mod (nrgs, NRGS) == 0)
+		    call realloc (rg, LEN_RG+2*(nrgs+NRGS), TY_STRUCT)
+		nrgs = nrgs + 1
+		RG_NRGS(rg) = nrgs
+		RG_X1(rg, nrgs) = max (rmin, rval1)
+		RG_X2(rg, nrgs) = min (rmax, rval2)
+		RG_NPTS(rg) = RG_NPTS(rg) +
+		    abs (RG_X1(rg, nrgs) - RG_X2(rg, nrgs)) + 1
+	    }
+	}
+
+	call sfree (sp)
+end
diff --git a/pkg/xtools/ranges/rgunion.x b/pkg/xtools/ranges/rgunion.x
new file mode 100644
index 00000000..5b9dfa6f
--- /dev/null
+++ b/pkg/xtools/ranges/rgunion.x
@@ -0,0 +1,48 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<pkg/rg.h>
+
+# RG_UNION -- Find the union of two sets of ranges.
+
+pointer procedure rg_union (rg1, rg2)
+
+pointer	rg1				# First set of ranges
+pointer	rg2				# Second set of ranges
+
+pointer	rg3				# Pointer to union
+
+int	i, j
+
+begin
+	# Error check the range pointers.
+
+	if ((rg1 == NULL) || (rg2 == NULL))
+	    call error (0, "Range descriptor(s) undefined")
+
+	# Allocate the range points array.
+
+	i = RG_NRGS(rg1) + RG_NRGS(rg2)
+	call malloc (rg3, LEN_RG + 2 * max (1, i), TY_STRUCT)
+
+	# Set the ranges.
+
+	RG_NRGS(rg3) = i
+	RG_NPTS(rg3) = RG_NPTS(rg1) + RG_NPTS(rg2)
+
+	j = 1
+	do i = 1, RG_NRGS(rg1) {
+	    RG_X1(rg3, j) = RG_X1(rg1, i)
+	    RG_X2(rg3, j) = RG_X2(rg1, i)
+	    j = j + 1
+	}
+	do i = 1, RG_NRGS(rg2) {
+	    RG_X1(rg3, j) = RG_X1(rg2, i)
+	    RG_X2(rg3, j) = RG_X2(rg2, i)
+	    j = j + 1
+	}
+
+	call rgorder (rg3)
+	call rgmerge (rg3)
+
+	return (rg3)
+end
diff --git a/pkg/xtools/ranges/rgunpack.gx b/pkg/xtools/ranges/rgunpack.gx
new file mode 100644
index 00000000..2b357ebb
--- /dev/null
+++ b/pkg/xtools/ranges/rgunpack.gx
@@ -0,0 +1,37 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<pkg/rg.h>
+
+# RG_UNPACK -- Unpack a packed array.
+#
+# There is no checking on the size of the arrays.  The points in the
+# unpacked array which are not covered by the packed array are left unchanged.
+# The packed and unpacked arrays should not be the same.
+
+procedure rg_unpack$t (rg, packed, unpacked)
+
+pointer	rg					# Ranges
+PIXEL	packed[ARB]				# Packed array
+PIXEL	unpacked[ARB]				# Unpacked array
+
+int	i, j, x1, x2, nx
+
+begin
+	if (rg == NULL)
+	    call error (0, "Range descriptor undefined")
+
+	j = 1
+	do i = 1, RG_NRGS(rg) {
+	    if (RG_X1(rg, i) < RG_X2(rg, i)) {
+		x1 = RG_X1(rg, i)
+		x2 = RG_X2(rg, i)
+	    } else {
+		x1 = RG_X2(rg, i)
+		x2 = RG_X1(rg, i)
+	    }
+
+	    nx = x2 - x1 + 1
+	    call amov$t (packed[j], unpacked[x1], nx)
+	    j = j + nx
+	}
+end
diff --git a/pkg/xtools/ranges/rgunpackd.x b/pkg/xtools/ranges/rgunpackd.x
new file mode 100644
index 00000000..2ce32fa2
--- /dev/null
+++ b/pkg/xtools/ranges/rgunpackd.x
@@ -0,0 +1,37 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<pkg/rg.h>
+
+# RG_UNPACK -- Unpack a packed array.
+#
+# There is no checking on the size of the arrays.  The points in the
+# unpacked array which are not covered by the packed array are left unchanged.
+# The packed and unpacked arrays should not be the same.
+
+procedure rg_unpackd (rg, packed, unpacked)
+
+pointer	rg					# Ranges
+double	packed[ARB]				# Packed array
+double	unpacked[ARB]				# Unpacked array
+
+int	i, j, x1, x2, nx
+
+begin
+	if (rg == NULL)
+	    call error (0, "Range descriptor undefined")
+
+	j = 1
+	do i = 1, RG_NRGS(rg) {
+	    if (RG_X1(rg, i) < RG_X2(rg, i)) {
+		x1 = RG_X1(rg, i)
+		x2 = RG_X2(rg, i)
+	    } else {
+		x1 = RG_X2(rg, i)
+		x2 = RG_X1(rg, i)
+	    }
+
+	    nx = x2 - x1 + 1
+	    call amovd (packed[j], unpacked[x1], nx)
+	    j = j + nx
+	}
+end
diff --git a/pkg/xtools/ranges/rgunpackr.x b/pkg/xtools/ranges/rgunpackr.x
new file mode 100644
index 00000000..6c96f5f8
--- /dev/null
+++ b/pkg/xtools/ranges/rgunpackr.x
@@ -0,0 +1,37 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<pkg/rg.h>
+
+# RG_UNPACK -- Unpack a packed array.
+#
+# There is no checking on the size of the arrays.  The points in the
+# unpacked array which are not covered by the packed array are left unchanged.
+# The packed and unpacked arrays should not be the same.
+
+procedure rg_unpackr (rg, packed, unpacked)
+
+pointer	rg					# Ranges
+real	packed[ARB]				# Packed array
+real	unpacked[ARB]				# Unpacked array
+
+int	i, j, x1, x2, nx
+
+begin
+	if (rg == NULL)
+	    call error (0, "Range descriptor undefined")
+
+	j = 1
+	do i = 1, RG_NRGS(rg) {
+	    if (RG_X1(rg, i) < RG_X2(rg, i)) {
+		x1 = RG_X1(rg, i)
+		x2 = RG_X2(rg, i)
+	    } else {
+		x1 = RG_X2(rg, i)
+		x2 = RG_X1(rg, i)
+	    }
+
+	    nx = x2 - x1 + 1
+	    call amovr (packed[j], unpacked[x1], nx)
+	    j = j + nx
+	}
+end
diff --git a/pkg/xtools/ranges/rgwindow.x b/pkg/xtools/ranges/rgwindow.x
new file mode 100644
index 00000000..fe495362
--- /dev/null
+++ b/pkg/xtools/ranges/rgwindow.x
@@ -0,0 +1,43 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<pkg/rg.h>
+
+# RG_WINDOW -- Intersect a set of ordered and merged ranges with a window.
+
+pointer procedure rg_window (rg, rmin, rmax)
+
+pointer	rg				# Ranges
+int	rmin, rmax			# Window
+
+pointer	rgout				# Pointer to windowed ranges
+
+int	i, j
+
+begin
+	if (rg == NULL)
+	    call error (0, "Range descriptor undefined")
+
+	# Allocate the range points array.
+
+	call malloc (rgout, LEN_RG + 2 * max (1, RG_NRGS(rg)), TY_STRUCT)
+
+	# Set the windowed ranges.
+
+	j = 0
+	do i = 1, RG_NRGS(rg) {
+	    if ((rmin <= RG_X2(rg, i)) && (rmax >= RG_X1(rg, i))) {
+		j = j + 1
+		RG_X1(rgout, j) = max (rmin, RG_X1(rg, i))
+		RG_X2(rgout, j) = min (rmax, RG_X2(rg, i))
+	    }
+	}
+
+	call realloc (rgout, LEN_RG + 2 * max (1, j), TY_STRUCT)
+	RG_NRGS(rgout) = j
+	RG_NPTS(rgout) = 0
+	do i = 1, RG_NRGS(rgout)
+	    RG_NPTS(rgout) = RG_NPTS(rgout) +
+		abs (RG_X1(rgout, i) - RG_X2(rgout, i)) + 1
+
+	return (rgout)
+end
diff --git a/pkg/xtools/ranges/rgwtbin.gx b/pkg/xtools/ranges/rgwtbin.gx
new file mode 100644
index 00000000..711dbf1e
--- /dev/null
+++ b/pkg/xtools/ranges/rgwtbin.gx
@@ -0,0 +1,112 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<pkg/rg.h>
+
+# RG_WTBIN -- Weighted average or median of data.
+#
+# The ranges are broken up into subranges of at most abs (nbin) points and a
+# minimum of max (3, (abs(nbin)+1)/2) (though always at least one bin).  The
+# subranges are weighted averaged if nbin > 1 and medianed if nbin < 1.
+# The output weights are the sum of the weights for each subrange.
+# The output array must be large enough to contain the desired points.
+# If the ranges are merged then the input and output arrays may be the same.
+
+procedure rg_wtbin$t (rg, nbin, in, wtin, nin, out, wtout, nout)
+
+pointer	rg				# Ranges
+int	nbin				# Maximum points in average or median
+PIXEL	in[nin]				# Input array
+PIXEL	wtin[nin]			# Input weights
+int	nin				# Number of input points
+PIXEL	out[ARB]			# Output array
+PIXEL	wtout[ARB]			# Output weights
+int	nout				# Number of output points
+
+int	i, j, k, l, n, npts, ntemp, nsample
+
+PIXEL	asum$t(), amed$t()
+
+errchk	rg_pack$t
+
+begin
+	# Check for a null set of ranges.
+
+	if (rg == NULL)
+	    call error (0, "Range descriptor undefined")
+
+	# If the bin size is exactly one then move the selected input points
+	# to the output array.
+
+	if (abs (nbin) < 2) {
+	    call rg_pack$t (rg, in, out)
+	    call rg_pack$t (rg, wtin, wtout)
+	    nout = RG_NPTS(rg)
+	    return
+	}
+
+	# Determine the subranges and take the median or average.
+
+	npts = abs (nbin)
+	ntemp = 0
+
+	do i = 1, RG_NRGS(rg) {
+	    nsample = 0
+	    if (RG_X1(rg, i) > RG_X2(rg, i)) {
+		j = min (nin, RG_X1(rg, i))
+		k = max (1, RG_X2(rg, i))
+		while (j >= k) {
+		    n = max (0, min (npts, j - k + 1))
+		    if (nsample > 0 && n < max (min (npts, 3), (npts+1)/2))
+			break
+		    k = k - n
+		    nsample = nsample + 1
+		    ntemp = ntemp + 1
+		    wtout[ntemp] = asum$t (wtin[k + 1], n)
+		    if (nbin > 0) {
+		        if (wtout[ntemp] != 0.) {
+			    out[ntemp] = 0.
+			    do l = k + 1, k + n
+			        out[ntemp] = out[ntemp] + in[l] * wtin[l]
+		            out[ntemp] = out[ntemp] / wtout[ntemp]
+			} else {
+			    out[ntemp] = 0.
+			    do l = k + 1, k + n
+			        out[ntemp] = out[ntemp] + in[l]
+		            out[ntemp] = out[ntemp] / n
+			}
+		    } else {
+			out[ntemp] = amed$t (in[k+1], n)
+		    }
+		}
+	    } else {
+		j = max (1, RG_X1(rg, i))
+		k = min (nin, RG_X2(rg, i))
+	        while (j <= k) {
+		    n = max (0, min (npts, k - j + 1))
+		    if (nsample > 0 && n < max (min (npts, 3), (npts+1)/2))
+			break
+		    nsample = nsample + 1
+		    ntemp = ntemp + 1
+		    wtout[ntemp] = asum$t (wtin[j], n)
+		    if (nbin > 0) {
+			if (wtout[ntemp] != 0.) {
+			    out[ntemp] = 0.
+			    do l = j, j + n - 1
+		                out[ntemp] = out[ntemp] + in[l] * wtin[l]
+			    out[ntemp] = out[ntemp] / wtout[ntemp]
+			} else {
+			    out[ntemp] = 0.
+			    do l = j, j + n - 1
+		                out[ntemp] = out[ntemp] + in[l]
+			    out[ntemp] = out[ntemp] / n
+			}
+		    } else {
+			out[ntemp] = amed$t (in[j], n)
+		    }
+		    j = j + n
+		}
+	    }
+	}
+
+	nout = ntemp
+end
diff --git a/pkg/xtools/ranges/rgwtbind.x b/pkg/xtools/ranges/rgwtbind.x
new file mode 100644
index 00000000..82adeba5
--- /dev/null
+++ b/pkg/xtools/ranges/rgwtbind.x
@@ -0,0 +1,112 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<pkg/rg.h>
+
+# RG_WTBIN -- Weighted average or median of data.
+#
+# The ranges are broken up into subranges of at most abs (nbin) points and a
+# minimum of max (3, (abs(nbin)+1)/2) (though always at least one bin).  The
+# subranges are weighted averaged if nbin > 1 and medianed if nbin < 1.
+# The output weights are the sum of the weights for each subrange.
+# The output array must be large enough to contain the desired points.
+# If the ranges are merged then the input and output arrays may be the same.
+
+procedure rg_wtbind (rg, nbin, in, wtin, nin, out, wtout, nout)
+
+pointer	rg				# Ranges
+int	nbin				# Maximum points in average or median
+double	in[nin]				# Input array
+double	wtin[nin]			# Input weights
+int	nin				# Number of input points
+double	out[ARB]			# Output array
+double	wtout[ARB]			# Output weights
+int	nout				# Number of output points
+
+int	i, j, k, l, n, npts, ntemp, nsample
+
+double	asumd(), amedd()
+
+errchk	rg_packd
+
+begin
+	# Check for a null set of ranges.
+
+	if (rg == NULL)
+	    call error (0, "Range descriptor undefined")
+
+	# If the bin size is exactly one then move the selected input points
+	# to the output array.
+
+	if (abs (nbin) < 2) {
+	    call rg_packd (rg, in, out)
+	    call rg_packd (rg, wtin, wtout)
+	    nout = RG_NPTS(rg)
+	    return
+	}
+
+	# Determine the subranges and take the median or average.
+
+	npts = abs (nbin)
+	ntemp = 0
+
+	do i = 1, RG_NRGS(rg) {
+	    nsample = 0
+	    if (RG_X1(rg, i) > RG_X2(rg, i)) {
+		j = min (nin, RG_X1(rg, i))
+		k = max (1, RG_X2(rg, i))
+		while (j >= k) {
+		    n = max (0, min (npts, j - k + 1))
+		    if (nsample > 0 && n < max (min (npts, 3), (npts+1)/2))
+			break
+		    k = k - n
+		    nsample = nsample + 1
+		    ntemp = ntemp + 1
+		    wtout[ntemp] = asumd (wtin[k + 1], n)
+		    if (nbin > 0) {
+		        if (wtout[ntemp] != 0.) {
+			    out[ntemp] = 0.
+			    do l = k + 1, k + n
+			        out[ntemp] = out[ntemp] + in[l] * wtin[l]
+		            out[ntemp] = out[ntemp] / wtout[ntemp]
+			} else {
+			    out[ntemp] = 0.
+			    do l = k + 1, k + n
+			        out[ntemp] = out[ntemp] + in[l]
+		            out[ntemp] = out[ntemp] / n
+			}
+		    } else {
+			out[ntemp] = amedd (in[k+1], n)
+		    }
+		}
+	    } else {
+		j = max (1, RG_X1(rg, i))
+		k = min (nin, RG_X2(rg, i))
+	        while (j <= k) {
+		    n = max (0, min (npts, k - j + 1))
+		    if (nsample > 0 && n < max (min (npts, 3), (npts+1)/2))
+			break
+		    nsample = nsample + 1
+		    ntemp = ntemp + 1
+		    wtout[ntemp] = asumd (wtin[j], n)
+		    if (nbin > 0) {
+			if (wtout[ntemp] != 0.) {
+			    out[ntemp] = 0.
+			    do l = j, j + n - 1
+		                out[ntemp] = out[ntemp] + in[l] * wtin[l]
+			    out[ntemp] = out[ntemp] / wtout[ntemp]
+			} else {
+			    out[ntemp] = 0.
+			    do l = j, j + n - 1
+		                out[ntemp] = out[ntemp] + in[l]
+			    out[ntemp] = out[ntemp] / n
+			}
+		    } else {
+			out[ntemp] = amedd (in[j], n)
+		    }
+		    j = j + n
+		}
+	    }
+	}
+
+	nout = ntemp
+end
diff --git a/pkg/xtools/ranges/rgwtbinr.x b/pkg/xtools/ranges/rgwtbinr.x
new file mode 100644
index 00000000..a4be8485
--- /dev/null
+++ b/pkg/xtools/ranges/rgwtbinr.x
@@ -0,0 +1,112 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<pkg/rg.h>
+
+# RG_WTBIN -- Weighted average or median of data.
+#
+# The ranges are broken up into subranges of at most abs (nbin) points and a
+# minimum of max (3, (abs(nbin)+1)/2) (though always at least one bin).  The
+# subranges are weighted averaged if nbin > 1 and medianed if nbin < 1.
+# The output weights are the sum of the weights for each subrange.
+# The output array must be large enough to contain the desired points.
+# If the ranges are merged then the input and output arrays may be the same.
+
+procedure rg_wtbinr (rg, nbin, in, wtin, nin, out, wtout, nout)
+
+pointer	rg				# Ranges
+int	nbin				# Maximum points in average or median
+real	in[nin]				# Input array
+real	wtin[nin]			# Input weights
+int	nin				# Number of input points
+real	out[ARB]			# Output array
+real	wtout[ARB]			# Output weights
+int	nout				# Number of output points
+
+int	i, j, k, l, n, npts, ntemp, nsample
+
+real	asumr(), amedr()
+
+errchk	rg_packr
+
+begin
+	# Check for a null set of ranges.
+
+	if (rg == NULL)
+	    call error (0, "Range descriptor undefined")
+
+	# If the bin size is exactly one then move the selected input points
+	# to the output array.
+
+	if (abs (nbin) < 2) {
+	    call rg_packr (rg, in, out)
+	    call rg_packr (rg, wtin, wtout)
+	    nout = RG_NPTS(rg)
+	    return
+	}
+
+	# Determine the subranges and take the median or average.
+
+	npts = abs (nbin)
+	ntemp = 0
+
+	do i = 1, RG_NRGS(rg) {
+	    nsample = 0
+	    if (RG_X1(rg, i) > RG_X2(rg, i)) {
+		j = min (nin, RG_X1(rg, i))
+		k = max (1, RG_X2(rg, i))
+		while (j >= k) {
+		    n = max (0, min (npts, j - k + 1))
+		    if (nsample > 0 && n < max (min (npts, 3), (npts+1)/2))
+			break
+		    k = k - n
+		    nsample = nsample + 1
+		    ntemp = ntemp + 1
+		    wtout[ntemp] = asumr (wtin[k + 1], n)
+		    if (nbin > 0) {
+		        if (wtout[ntemp] != 0.) {
+			    out[ntemp] = 0.
+			    do l = k + 1, k + n
+			        out[ntemp] = out[ntemp] + in[l] * wtin[l]
+		            out[ntemp] = out[ntemp] / wtout[ntemp]
+			} else {
+			    out[ntemp] = 0.
+			    do l = k + 1, k + n
+			        out[ntemp] = out[ntemp] + in[l]
+		            out[ntemp] = out[ntemp] / n
+			}
+		    } else {
+			out[ntemp] = amedr (in[k+1], n)
+		    }
+		}
+	    } else {
+		j = max (1, RG_X1(rg, i))
+		k = min (nin, RG_X2(rg, i))
+	        while (j <= k) {
+		    n = max (0, min (npts, k - j + 1))
+		    if (nsample > 0 && n < max (min (npts, 3), (npts+1)/2))
+			break
+		    nsample = nsample + 1
+		    ntemp = ntemp + 1
+		    wtout[ntemp] = asumr (wtin[j], n)
+		    if (nbin > 0) {
+			if (wtout[ntemp] != 0.) {
+			    out[ntemp] = 0.
+			    do l = j, j + n - 1
+		                out[ntemp] = out[ntemp] + in[l] * wtin[l]
+			    out[ntemp] = out[ntemp] / wtout[ntemp]
+			} else {
+			    out[ntemp] = 0.
+			    do l = j, j + n - 1
+		                out[ntemp] = out[ntemp] + in[l]
+			    out[ntemp] = out[ntemp] / n
+			}
+		    } else {
+			out[ntemp] = amedr (in[j], n)
+		    }
+		    j = j + n
+		}
+	    }
+	}
+
+	nout = ntemp
+end
diff --git a/pkg/xtools/ranges/rgxranges.gx b/pkg/xtools/ranges/rgxranges.gx
new file mode 100644
index 00000000..7a779925
--- /dev/null
+++ b/pkg/xtools/ranges/rgxranges.gx
@@ -0,0 +1,162 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<error.h>
+include	<ctype.h>
+include	<pkg/rg.h>
+
+define	NRGS	10		# Allocation size
+
+# RG_XRANGES -- Parse a range string corrsponding to a real set of values.
+# Return a pointer to the ranges.
+
+pointer procedure rg_xranges$t (rstr, rvals, npts)
+
+char	rstr[ARB]		# Range string
+PIXEL	rvals[npts]		# Range values (sorted)
+int	npts			# Number of range values
+pointer	rg			# Range pointer
+
+int	i, fd, strlen(), open(), getline()
+pointer	sp, str, ptr
+errchk	open, rg_xadd$t
+
+begin
+	# Check for valid arguments
+	if (npts < 1)
+	    call error (0, "No data points for range determination")
+
+	call smark (sp)
+	call salloc (str, max (strlen (rstr), SZ_LINE), TY_CHAR)
+	call calloc (rg, LEN_RG, TY_STRUCT)
+
+	i = 1
+	while (rstr[i] != EOS) {
+
+	    # Find beginning and end of a range and copy it to the work string
+	    while (IS_WHITE(rstr[i]) || rstr[i]==',' || rstr[i]=='\n')
+	        i = i + 1
+	    if (rstr[i] == EOS)
+		break
+
+	    ptr = str
+	    while (!(IS_WHITE(rstr[i]) || rstr[i]==',' || rstr[i]=='\n' ||
+		rstr[i]==EOS)) {
+		Memc[ptr] = rstr[i]
+	        i = i + 1
+		ptr = ptr + 1
+	    }
+	    Memc[ptr] = EOS
+
+	    # Add range(s)
+	    iferr {
+		if (Memc[str] == '@') {
+		    fd = open (Memc[str+1], READ_ONLY, TEXT_FILE)
+		    while (getline (fd, Memc[str]) != EOF) {
+			iferr (call rg_xadd$t (rg, Memc[str], rvals, npts))
+			    call erract (EA_WARN)
+		    }
+		    call close (fd)
+		} else
+		    call rg_xadd$t (rg, Memc[str], rvals, npts)
+	    } then
+		call erract (EA_WARN)
+	}
+
+	call sfree (sp)
+	return (rg)
+end
+
+
+# RG_XADD -- Add a range
+
+procedure rg_xadd$t (rg, rstr, rvals, npts)
+
+pointer	rg			# Range descriptor
+char	rstr[ARB]		# Range string
+PIXEL	rvals[npts]		# Range values (sorted)
+int	npts			# Number of range values
+
+int	i, j, k, nrgs, strlen(), cto$t()
+PIXEL	rval1, rval2, a1, b1, a2, b2
+pointer	sp, str, ptr
+
+begin
+	call smark (sp)
+	call salloc (str, strlen (rstr), TY_CHAR)
+
+	i = 1
+	while (rstr[i] != EOS) {
+
+	    # Find beginning and end of a range and copy it to the work string
+	    while (IS_WHITE(rstr[i]) || rstr[i]==',' || rstr[i]=='\n')
+	        i = i + 1
+	    if (rstr[i] == EOS)
+		break
+
+	    ptr = str
+	    while (!(IS_WHITE(rstr[i]) || rstr[i]==',' || rstr[i]=='\n' ||
+		rstr[i]==EOS)) {
+		if (rstr[i] == ':')
+		    Memc[ptr] = ' '
+		else
+		    Memc[ptr] = rstr[i]
+	        i = i + 1
+		ptr = ptr + 1
+	    }
+	    Memc[ptr] = EOS
+
+	    # Parse range
+	    if (Memc[str] == '@')
+		call error (1, "Cannot nest @files")
+	    else if (Memc[str] == '*') {
+		rval1 = rvals[1]
+		rval2 = rvals[npts]
+	    } else {
+		# Get range
+		j = 1
+		if (cto$t (Memc[str], j, rval1) == 0)
+		    call error (1, "Range syntax error")
+		rval2 = rval1
+		if (cto$t (Memc[str], j, rval2) == 0)
+		    ;
+	    }
+
+	    # Check limits and find indices into rval array
+	    a1 = min (rval1, rval2)
+	    b1 = max (rval1, rval2)
+	    a2 = min (rvals[1], rvals[npts])
+	    b2 = max (rvals[1], rvals[npts])
+	    if ((b1 >= a2) && (a1 <= b2)) {
+		a1 = max (a2, min (b2, a1))
+		b1 = max (a2, min (b2, b1))
+		if (rvals[1] <= rvals[npts]) {
+		    for (k = 1; (k <= npts) && (rvals[k] < a1); k = k + 1)
+			;
+		    for (j = k; (j <= npts) && (rvals[j] <= b1); j = j + 1)
+			;
+		    j = j - 1
+		} else {
+		    for (k = 1; (k <= npts) && (rvals[k] > b1); k = k + 1)
+			;
+		    for (j = k; (j <= npts) && (rvals[j] >= a1); j = j + 1)
+			;
+		    j = j - 1
+		}
+
+		# Add range
+		if (k <= j) {
+		    nrgs = RG_NRGS(rg)
+		    if (mod (nrgs, NRGS) == 0)
+			call realloc (rg, LEN_RG+2*(nrgs+NRGS), TY_STRUCT)
+		    nrgs = nrgs + 1
+		    RG_NRGS(rg) = nrgs
+		    RG_X1(rg, nrgs) = k
+		    RG_X2(rg, nrgs) = j
+		    RG_NPTS(rg) = RG_NPTS(rg) +
+			RG_X1(rg, nrgs) - RG_X2(rg, nrgs) + 1
+		}
+	    }
+	}
+
+	call sfree (sp)
+end
diff --git a/pkg/xtools/ranges/rgxranges1.gx b/pkg/xtools/ranges/rgxranges1.gx
new file mode 100644
index 00000000..b019e47c
--- /dev/null
+++ b/pkg/xtools/ranges/rgxranges1.gx
@@ -0,0 +1,146 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<ctype.h>
+include	<pkg/rg.h>
+
+# RG_XRANGES -- Parse a range string corrsponding to a real set of values.
+# Return a pointer to the ranges.
+
+pointer procedure rg_xranges$t (rstr, rvals, npts)
+
+char	rstr[ARB]		# Range string
+PIXEL	rvals[npts]		# Range values (sorted)
+int	npts			# Number of range values
+
+pointer	rg
+int	i, j, k, nrgs
+PIXEL	rval1, rval2, a, b
+
+int	cto$t()
+
+begin
+	# Check for valid arguments.
+
+	if (npts < 1)
+	    call error (0, "No data points for range determination")
+
+	# Check for a valid string and determine the number of ranges.
+
+	i = 1
+	nrgs = 0
+
+	while (rstr[i] != EOS) {
+
+	    # Skip delimiters
+	    while (IS_WHITE(rstr[i]) || (rstr[i] == ',') || (rstr[i]=='\n'))
+	        i = i + 1
+
+	    # Check for EOS.
+
+	    if (rstr[i] == EOS)
+		break
+
+	    # First character must be a *, -, ., or digit.
+
+	    if ((rstr[i] == '*') || (rstr[i] == '-') || (rstr[i] == '.') ||
+		    IS_DIGIT(rstr[i])) {
+		i = i + 1
+	        nrgs = nrgs + 1
+
+ 		# Remaining characters must be :, -, ., E, D, e, d, or digits.
+		# Replace : with ! to avoid sexigesimal interpretation.
+
+	        while ((rstr[i]==':') || (rstr[i]=='-') || (rstr[i]=='.') ||
+ 		       (rstr[i]=='E') || (rstr[i]=='D') ||
+ 		       (rstr[i]=='e') || (rstr[i]=='d') ||
+			IS_DIGIT(rstr[i])) {
+		    if (rstr[i] == ':')
+			rstr[i] = '!'
+		    i = i + 1
+		}
+	    } else
+		call error (0, "Syntax error in range string")
+	}
+
+	# Allocate memory for the ranges.
+
+	call malloc (rg, LEN_RG + 2 * max (1, nrgs), TY_STRUCT)
+
+	# Rescan the string and set the ranges.
+
+	i = 1
+	nrgs = 0
+
+	while (rstr[i] != EOS) {
+
+	    # Skip delimiters.
+	    while (IS_WHITE(rstr[i]) || (rstr[i]==',') || (rstr[i]=='\n'))
+	        i = i + 1
+
+	    # Check for EOS.
+
+	    if (rstr[i] == EOS)
+		break
+
+	    # If first character is * then set range to full range.
+	    # Otherwise parse the range.
+
+	    if (rstr[i] == '*') {
+		i = i + 1
+		rval1 = rvals[1]
+		rval2 = rvals[npts]
+
+	    } else {
+		# First digit is starting value.
+		if (cto$t (rstr, i, rval1) == 0) {
+		    nrgs = 0
+		    break
+		}
+		rval2 = rval1
+
+		# Check for an ending value for the range and restore ':'.
+		if (rstr[i] == '!') {
+		    rstr[i] = ':'
+		    i = i + 1
+		    if (cto$t (rstr, i, rval2) == 0) {
+			nrgs = 0
+			break
+		    }
+		}
+	    }
+
+	    # Check limits.
+
+	    a = min (rval1, rval2)
+	    b = max (rval1, rval2)
+	    if ((b >= rvals[1]) && (a <= rvals[npts])) {
+		rval1 = max (rvals[1], min (rvals[npts], rval1))
+		rval2 = max (rvals[1], min (rvals[npts], rval2))
+	    	a = min (rval1, rval2)
+	    	b = max (rval1, rval2)
+		for (k = 1; (k <= npts) && (rvals[k] < a); k = k + 1)
+		    ;
+		for (j = k; (j <= npts) && (rvals[j] <= b); j = j + 1)
+		    ;
+		j = j - 1
+		if (k <= j) {
+		    nrgs = nrgs + 1
+		    if (rval1 <= rval2) {
+		    	RG_X1(rg, nrgs) = k
+		    	RG_X2(rg, nrgs) = j
+		    } else {
+		    	RG_X1(rg, nrgs) = j
+		    	RG_X2(rg, nrgs) = k
+		    }
+		}
+	    }
+	}
+
+	RG_NRGS(rg) = nrgs
+	RG_NPTS(rg) = 0
+	do i = 1, RG_NRGS(rg)
+	    RG_NPTS(rg) = RG_NPTS(rg) +
+		abs (RG_X1(rg, i) - RG_X2(rg, i)) + 1
+
+	return (rg)
+end
diff --git a/pkg/xtools/ranges/rgxrangesd.x b/pkg/xtools/ranges/rgxrangesd.x
new file mode 100644
index 00000000..f9de6c32
--- /dev/null
+++ b/pkg/xtools/ranges/rgxrangesd.x
@@ -0,0 +1,162 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<error.h>
+include	<ctype.h>
+include	<pkg/rg.h>
+
+define	NRGS	10		# Allocation size
+
+# RG_XRANGES -- Parse a range string corrsponding to a real set of values.
+# Return a pointer to the ranges.
+
+pointer procedure rg_xrangesd (rstr, rvals, npts)
+
+char	rstr[ARB]		# Range string
+double	rvals[npts]		# Range values (sorted)
+int	npts			# Number of range values
+pointer	rg			# Range pointer
+
+int	i, fd, strlen(), open(), getline()
+pointer	sp, str, ptr
+errchk	open, rg_xaddd
+
+begin
+	# Check for valid arguments
+	if (npts < 1)
+	    call error (0, "No data points for range determination")
+
+	call smark (sp)
+	call salloc (str, max (strlen (rstr), SZ_LINE), TY_CHAR)
+	call calloc (rg, LEN_RG, TY_STRUCT)
+
+	i = 1
+	while (rstr[i] != EOS) {
+
+	    # Find beginning and end of a range and copy it to the work string
+	    while (IS_WHITE(rstr[i]) || rstr[i]==',' || rstr[i]=='\n')
+	        i = i + 1
+	    if (rstr[i] == EOS)
+		break
+
+	    ptr = str
+	    while (!(IS_WHITE(rstr[i]) || rstr[i]==',' || rstr[i]=='\n' ||
+		rstr[i]==EOS)) {
+		Memc[ptr] = rstr[i]
+	        i = i + 1
+		ptr = ptr + 1
+	    }
+	    Memc[ptr] = EOS
+
+	    # Add range(s)
+	    iferr {
+		if (Memc[str] == '@') {
+		    fd = open (Memc[str+1], READ_ONLY, TEXT_FILE)
+		    while (getline (fd, Memc[str]) != EOF) {
+			iferr (call rg_xaddd (rg, Memc[str], rvals, npts))
+			    call erract (EA_WARN)
+		    }
+		    call close (fd)
+		} else
+		    call rg_xaddd (rg, Memc[str], rvals, npts)
+	    } then
+		call erract (EA_WARN)
+	}
+
+	call sfree (sp)
+	return (rg)
+end
+
+
+# RG_XADD -- Add a range
+
+procedure rg_xaddd (rg, rstr, rvals, npts)
+
+pointer	rg			# Range descriptor
+char	rstr[ARB]		# Range string
+double	rvals[npts]		# Range values (sorted)
+int	npts			# Number of range values
+
+int	i, j, k, nrgs, strlen(), ctod()
+double	rval1, rval2, a1, b1, a2, b2
+pointer	sp, str, ptr
+
+begin
+	call smark (sp)
+	call salloc (str, strlen (rstr), TY_CHAR)
+
+	i = 1
+	while (rstr[i] != EOS) {
+
+	    # Find beginning and end of a range and copy it to the work string
+	    while (IS_WHITE(rstr[i]) || rstr[i]==',' || rstr[i]=='\n')
+	        i = i + 1
+	    if (rstr[i] == EOS)
+		break
+
+	    ptr = str
+	    while (!(IS_WHITE(rstr[i]) || rstr[i]==',' || rstr[i]=='\n' ||
+		rstr[i]==EOS)) {
+		if (rstr[i] == ':')
+		    Memc[ptr] = ' '
+		else
+		    Memc[ptr] = rstr[i]
+	        i = i + 1
+		ptr = ptr + 1
+	    }
+	    Memc[ptr] = EOS
+
+	    # Parse range
+	    if (Memc[str] == '@')
+		call error (1, "Cannot nest @files")
+	    else if (Memc[str] == '*') {
+		rval1 = rvals[1]
+		rval2 = rvals[npts]
+	    } else {
+		# Get range
+		j = 1
+		if (ctod (Memc[str], j, rval1) == 0)
+		    call error (1, "Range syntax error")
+		rval2 = rval1
+		if (ctod (Memc[str], j, rval2) == 0)
+		    ;
+	    }
+
+	    # Check limits and find indices into rval array
+	    a1 = min (rval1, rval2)
+	    b1 = max (rval1, rval2)
+	    a2 = min (rvals[1], rvals[npts])
+	    b2 = max (rvals[1], rvals[npts])
+	    if ((b1 >= a2) && (a1 <= b2)) {
+		a1 = max (a2, min (b2, a1))
+		b1 = max (a2, min (b2, b1))
+		if (rvals[1] <= rvals[npts]) {
+		    for (k = 1; (k <= npts) && (rvals[k] < a1); k = k + 1)
+			;
+		    for (j = k; (j <= npts) && (rvals[j] <= b1); j = j + 1)
+			;
+		    j = j - 1
+		} else {
+		    for (k = 1; (k <= npts) && (rvals[k] > b1); k = k + 1)
+			;
+		    for (j = k; (j <= npts) && (rvals[j] >= a1); j = j + 1)
+			;
+		    j = j - 1
+		}
+
+		# Add range
+		if (k <= j) {
+		    nrgs = RG_NRGS(rg)
+		    if (mod (nrgs, NRGS) == 0)
+			call realloc (rg, LEN_RG+2*(nrgs+NRGS), TY_STRUCT)
+		    nrgs = nrgs + 1
+		    RG_NRGS(rg) = nrgs
+		    RG_X1(rg, nrgs) = k
+		    RG_X2(rg, nrgs) = j
+		    RG_NPTS(rg) = RG_NPTS(rg) +
+			RG_X1(rg, nrgs) - RG_X2(rg, nrgs) + 1
+		}
+	    }
+	}
+
+	call sfree (sp)
+end
diff --git a/pkg/xtools/ranges/rgxrangesr.x b/pkg/xtools/ranges/rgxrangesr.x
new file mode 100644
index 00000000..425abf04
--- /dev/null
+++ b/pkg/xtools/ranges/rgxrangesr.x
@@ -0,0 +1,162 @@
+# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
+
+include	<error.h>
+include	<ctype.h>
+include	<pkg/rg.h>
+
+define	NRGS	10		# Allocation size
+
+# RG_XRANGES -- Parse a range string corrsponding to a real set of values.
+# Return a pointer to the ranges.
+
+pointer procedure rg_xrangesr (rstr, rvals, npts)
+
+char	rstr[ARB]		# Range string
+real	rvals[npts]		# Range values (sorted)
+int	npts			# Number of range values
+pointer	rg			# Range pointer
+
+int	i, fd, strlen(), open(), getline()
+pointer	sp, str, ptr
+errchk	open, rg_xaddr
+
+begin
+	# Check for valid arguments
+	if (npts < 1)
+	    call error (0, "No data points for range determination")
+
+	call smark (sp)
+	call salloc (str, max (strlen (rstr), SZ_LINE), TY_CHAR)
+	call calloc (rg, LEN_RG, TY_STRUCT)
+
+	i = 1
+	while (rstr[i] != EOS) {
+
+	    # Find beginning and end of a range and copy it to the work string
+	    while (IS_WHITE(rstr[i]) || rstr[i]==',' || rstr[i]=='\n')
+	        i = i + 1
+	    if (rstr[i] == EOS)
+		break
+
+	    ptr = str
+	    while (!(IS_WHITE(rstr[i]) || rstr[i]==',' || rstr[i]=='\n' ||
+		rstr[i]==EOS)) {
+		Memc[ptr] = rstr[i]
+	        i = i + 1
+		ptr = ptr + 1
+	    }
+	    Memc[ptr] = EOS
+
+	    # Add range(s)
+	    iferr {
+		if (Memc[str] == '@') {
+		    fd = open (Memc[str+1], READ_ONLY, TEXT_FILE)
+		    while (getline (fd, Memc[str]) != EOF) {
+			iferr (call rg_xaddr (rg, Memc[str], rvals, npts))
+			    call erract (EA_WARN)
+		    }
+		    call close (fd)
+		} else
+		    call rg_xaddr (rg, Memc[str], rvals, npts)
+	    } then
+		call erract (EA_WARN)
+	}
+
+	call sfree (sp)
+	return (rg)
+end
+
+
+# RG_XADD -- Add a range
+
+procedure rg_xaddr (rg, rstr, rvals, npts)
+
+pointer	rg			# Range descriptor
+char	rstr[ARB]		# Range string
+real	rvals[npts]		# Range values (sorted)
+int	npts			# Number of range values
+
+int	i, j, k, nrgs, strlen(), ctor()
+real	rval1, rval2, a1, b1, a2, b2
+pointer	sp, str, ptr
+
+begin
+	call smark (sp)
+	call salloc (str, strlen (rstr), TY_CHAR)
+
+	i = 1
+	while (rstr[i] != EOS) {
+
+	    # Find beginning and end of a range and copy it to the work string
+	    while (IS_WHITE(rstr[i]) || rstr[i]==',' || rstr[i]=='\n')
+	        i = i + 1
+	    if (rstr[i] == EOS)
+		break
+
+	    ptr = str
+	    while (!(IS_WHITE(rstr[i]) || rstr[i]==',' || rstr[i]=='\n' ||
+		rstr[i]==EOS)) {
+		if (rstr[i] == ':')
+		    Memc[ptr] = ' '
+		else
+		    Memc[ptr] = rstr[i]
+	        i = i + 1
+		ptr = ptr + 1
+	    }
+	    Memc[ptr] = EOS
+
+	    # Parse range
+	    if (Memc[str] == '@')
+		call error (1, "Cannot nest @files")
+	    else if (Memc[str] == '*') {
+		rval1 = rvals[1]
+		rval2 = rvals[npts]
+	    } else {
+		# Get range
+		j = 1
+		if (ctor (Memc[str], j, rval1) == 0)
+		    call error (1, "Range syntax error")
+		rval2 = rval1
+		if (ctor (Memc[str], j, rval2) == 0)
+		    ;
+	    }
+
+	    # Check limits and find indices into rval array
+	    a1 = min (rval1, rval2)
+	    b1 = max (rval1, rval2)
+	    a2 = min (rvals[1], rvals[npts])
+	    b2 = max (rvals[1], rvals[npts])
+	    if ((b1 >= a2) && (a1 <= b2)) {
+		a1 = max (a2, min (b2, a1))
+		b1 = max (a2, min (b2, b1))
+		if (rvals[1] <= rvals[npts]) {
+		    for (k = 1; (k <= npts) && (rvals[k] < a1); k = k + 1)
+			;
+		    for (j = k; (j <= npts) && (rvals[j] <= b1); j = j + 1)
+			;
+		    j = j - 1
+		} else {
+		    for (k = 1; (k <= npts) && (rvals[k] > b1); k = k + 1)
+			;
+		    for (j = k; (j <= npts) && (rvals[j] >= a1); j = j + 1)
+			;
+		    j = j - 1
+		}
+
+		# Add range
+		if (k <= j) {
+		    nrgs = RG_NRGS(rg)
+		    if (mod (nrgs, NRGS) == 0)
+			call realloc (rg, LEN_RG+2*(nrgs+NRGS), TY_STRUCT)
+		    nrgs = nrgs + 1
+		    RG_NRGS(rg) = nrgs
+		    RG_X1(rg, nrgs) = k
+		    RG_X2(rg, nrgs) = j
+		    RG_NPTS(rg) = RG_NPTS(rg) +
+			RG_X1(rg, nrgs) - RG_X2(rg, nrgs) + 1
+		}
+	    }
+	}
+
+	call sfree (sp)
+end