Initial commit

47 files changed, 5434 insertions, 0 deletions

diff --git a/noao/onedspec/irsiids/addsets.par b/noao/onedspec/irsiids/addsets.par
new file mode 100644
index 00000000..1ad84e8b
--- /dev/null
+++ b/noao/onedspec/irsiids/addsets.par
@@ -0,0 +1,8 @@
+# ADDSETS parameter file
+
+input,s,a,,,,Input image root file name
+records,s,a,,,,Range of spectral records
+output,s,a,,,,Output file root name for new spectra
+start_rec,i,a,1,0,9999,Next starting spectral record
+subset,i,h,2,1,,Number of spectra to add together within string
+weighting,b,h,yes,,,Apply integration time weighting to calibrated data
diff --git a/noao/onedspec/irsiids/batchred.cl b/noao/onedspec/irsiids/batchred.cl
new file mode 100644
index 00000000..1bd5a3a3
--- /dev/null
+++ b/noao/onedspec/irsiids/batchred.cl
@@ -0,0 +1,168 @@
+#{ BATCHRED -- Script file to generate another script file
+#  which runs several ONEDSPEC tasks in an automated fashion.
+#
+#  Currently the following procedures are automated:
+#    1. STANDARD
+#    2. SENSFUNC
+#    3. BSWITCH
+#    4. CALIBRATE
+#    5. ADDSETS
+#
+
+{
+# Say hello to the guy on the other side of the screen and check batch file.
+print ("\n----B A T C H   I I D S / I R S   F I L E   G E N E R A T O R----\n")
+
+s2 = "process.cl"			# Batch file to be created.
+if (access (s2)) {
+    print ("A batch file already exists - ")
+    if (query)
+        delete (s2, verify=no)
+}
+
+# Initialize
+rt = input				# Root name for spectra
+ot = output				# Output root name
+ttl = ">>&'" + ttylog + "')\n"		# Log file for tty output
+
+out = rt
+s1 = ""
+st = ""
+sns = ""
+stat = ""
+print ("i = ", start_rec, >>s2)
+
+if (standard) {				# STANDARD?
+    print ("\n#---STANDARD---\n")
+    print ("\n#---STANDARD---\n", >>s2)
+
+    st = std				# STD file
+    if (access (st)) {
+	print (st, " - already exists")
+	if (query)
+            delete (st, verify=no)
+    }
+
+    # Loop over all stars
+    b1 = yes
+    while (b1) {
+        # Check that the last entry was different - otherwise end input
+	records = ""
+    	s3 = records
+    	if (s3 == "")
+    	    b1 = no
+        else {
+    	    print ("standard (input='",rt,"',output='",st,"',",>>s2)
+    	    print ("\trecords='",s3,"',",>>s2)
+    	    print ("\tstar_name='",star_name,"',beam_switch=yes,",>>s2)
+	    print ("\tsamestar=yes,apertures='',bandwidth=INDEF,",>>s2)
+	    print ("\tbandsep=INDEF,interact=no,",ttl,>>s2)
+        }
+    }
+
+    print ("")
+}
+
+
+if (sensfunc) {				# SENSFUNC?
+    print ("\n#---SENSFUNC---\n")
+    print ("\n#---SENSFUNC---\n", >>s2)
+
+    if (st == "")
+	st = std			# STD file
+    sns = sensitivity			# Sensitivity image
+    stat = stats			# Statistics file
+
+    print ("\nsensfunc (standards='",st,"',sensitivity='",sns,"',",>>s2)
+    print ("\tlogfile='",stat,"',apertures='',ignoreaps=no,",>>s2)
+    print ("\tfunction='",function,"',order=",order,",",>>s2)
+    print ("\tinteract=no,",ttl,>>s2)
+
+    print ("")
+}
+
+
+if (bswitch) {				# BSWITCH?
+    print ("\n#---BSWITCH---\n")
+    print ("\n#---BSWITCH---\n", >>s2)
+
+    # Save starting output record number
+    in = out
+    out = "b" // ot
+    wt = weight				# Weighting?
+    if (stat == "")
+	stat = stats			# Statistics file
+
+    # Accumulate records
+    print ("next_rec = i", >>s2)
+    b1 = yes
+    while (b1) {
+	records = ""
+        s3 = records
+        if (s3 == "")
+    	    b1 = no
+        else {
+	    print ("j = next_rec", >> s2)
+    	    print ("bswitch (input='",in,"',output='",out,"',",>>s2)
+    	    print ("\trecords='",s3,"',stats='",stat,"',",>>s2)
+    	    print ("\tweighting=",wt,",subset=",subset,",",>>s2)
+	    print ("\tstart_rec=j,",>>s2)
+    	    print ("\twave1=",wave1,",wave2=",wave2,",",ttl,>>s2)
+        }
+    }
+
+    # Output records
+    print ("j = next_rec", >>s2)
+    s1 = "str (i) // '-' // str(j-1)"
+    print ("s1 = ", s1, >>s2)
+
+    print ("")
+}
+
+if (calibrate) {			# CALIBRATE?
+    print ("\n#---CALIBRATE---\n")
+    print ("\n#---CALIBRATE---\n", >>s2)
+
+    in = out
+    out = "c" // ot
+    if (sns == "")
+        sns = sensitivity		# Sensivity file name
+
+    if (s1 == "") {
+	records = ""
+	s1 = records
+        print ("s1 = '", s1, "'", >>s2)
+    }
+
+    print ("calibrate (input='",in,"',output='",out,"',records=s1,",>>s2)
+    print ("\tignoreaps=no,",>>s2)
+    print ("\textinct=no,flux=yes,",>>s2)
+    print ("\tsensitivity='",sns,"',fnu=",fnu,",",ttl,>>s2)
+
+    print ("")
+}
+
+if (addsets) {				# ADDSETS?
+    print ("\n#---ADDSETS---\n")
+    print ("\n#---ADDSETS---\n", >>s2)
+
+    in = out
+    out = "a" // ot
+    if (s1 == "") {
+	records = ""
+	s1 = records
+        print ("s1 = '", s1, "'", >>s2)
+    }
+
+    print ("addsets (input='",in,"',output='",out,"',records=s1,",>>s2)
+    print ("\tstart_rec=i,subset=2,",ttl,>>s2)
+}
+
+# All done with generator.  Ask whether to execute it.
+print ("File generation complete - filename=",s2)
+if (proceed == no)			# Execute batch file?
+    bye
+}
+
+# Execute generated batch file
+process &
diff --git a/noao/onedspec/irsiids/batchred.par b/noao/onedspec/irsiids/batchred.par
new file mode 100644
index 00000000..dbb61b8c
--- /dev/null
+++ b/noao/onedspec/irsiids/batchred.par
@@ -0,0 +1,38 @@
+# BATCHRED -- Parameter file for batch reduction prep task
+
+input,s,a,,,,Input root name for spectra
+output,s,a,,,,Output root name for spectra
+start_rec,i,a,,0,9999,Starting record for output spectra
+ttylog,s,a,"ttylog",,,File name to contain a log of terminal output
+standard,b,a,yes,,,Generate commands for STANDARD
+sensfunc,b,a,yes,,,Generate commands for SENSFUNC
+bswitch,b,a,yes,,,Generate commands for BSWITCH
+calibrate,b,a,yes,,,Generate commands for CALIBRATE
+addsets,b,a,yes,,,Generate commands for ADDSETS
+
+std,s,a,"std",,,STANDARD and SENSFUNC standard star file
+star_name,s,q,,,,STANDARD star name
+stats,s,a,"stats",,,SENSFUNC and BSWITCH statistics file
+sensitivity,s,a,sens,,,SENSFUNC and CALIBRATE sensitivity spectra
+weight,b,a,no,,,BSWITCH weighted averages?
+function,s,h,"chebyshev",,,SENSFUNC fitting function
+order,i,h,7,1,,SENSFUNC fitting order
+
+records,s,q,,,,Record string to process
+proceed,b,q,yes,,,Begin batch processing
+query,b,q,no,,,Delete files(s)?
+
+fnu,b,h,no
+wave1,r,h,0.0
+wave2,r,h,0.0
+subset,i,h,32767
+
+rt,s,h
+ot,s,h
+in,s,h
+out,s,h
+stat,s,h
+sns,s,h
+st,s,h
+wt,b,h
+ttl,s,h
diff --git a/noao/onedspec/irsiids/bplot.cl b/noao/onedspec/irsiids/bplot.cl
new file mode 100644
index 00000000..53ce4cfc
--- /dev/null
+++ b/noao/onedspec/irsiids/bplot.cl
@@ -0,0 +1,35 @@
+# BPLOT -- Batch plotting of spectra with SPLOT
+
+procedure bplot (images, records)
+
+string	images				{prompt="List of images to plot"}
+string	records = ""			{prompt="List of records to plot"}
+string	graphics = "stdgraph"		{prompt="Graphics output device"}
+string	cursor   = "onedspec$gcurval.dat"	{prompt="Cursor file(s)\n"}
+
+struct	*ilist, *clist
+
+begin
+	int	line, ap
+	file	ifile, cfile, cur, image
+
+	ifile = mktemp ("bplot")
+	cfile = mktemp ("bplot")
+
+	names (images, records, >& ifile)
+	files (cursor, > cfile) 
+	cur = ""
+
+	ilist = ifile; clist = cfile
+	while (fscan (ilist, image) != EOF) {
+	    if ((cursor != "") && (fscan (clist, cur) == EOF)) {
+		clist = cfile
+		line = fscan (clist, cur)
+	    }
+	    splot (image, graphics=graphics, cursor=cur)
+	}
+	clist = ""; ilist = ""
+
+	delete (ifile, verify=no)
+	delete (cfile, verify=no)
+end
diff --git a/noao/onedspec/irsiids/bswitch.par b/noao/onedspec/irsiids/bswitch.par
new file mode 100644
index 00000000..3a8b7e8f
--- /dev/null
+++ b/noao/onedspec/irsiids/bswitch.par
@@ -0,0 +1,15 @@
+# BSWITCH parameter file
+
+input,s,a,,,,Input spectra file root name
+records,s,a,,,,Ranges of spectral records
+output,s,a,,,,Output file root name for new spectra
+start_rec,i,a,1,0,9999,Next starting spectral record
+stats,s,a,"stats",,,File to contain statistics summary
+ids_mode,b,h,yes,,,Are data in quadruples
+extinct,b,h,yes,,,Apply de-extinction correction
+weighting,b,h,no,,,Apply statistical weights during addition
+subset,i,h,32767,1,,Generate sums at subset intervals
+wave1,r,h,0.0,,,Starting wavelength to accumulate stats
+wave2,r,h,0.0,,,Ending wavelength to accumulate stats
+observatory,s,h,"kpno",,,Observatory of data
+extinction,s,h,)_.extinction,,,Extinction file
diff --git a/noao/onedspec/irsiids/coefs.par b/noao/onedspec/irsiids/coefs.par
new file mode 100644
index 00000000..dbabec65
--- /dev/null
+++ b/noao/onedspec/irsiids/coefs.par
@@ -0,0 +1,3 @@
+input,s,a,,,,Input image root file name
+records,s,a,,,,Range of spectral records
+database,s,h,"database",,,IDENTIFY database
diff --git a/noao/onedspec/irsiids/coincor.par b/noao/onedspec/irsiids/coincor.par
new file mode 100644
index 00000000..a0f2e0bb
--- /dev/null
+++ b/noao/onedspec/irsiids/coincor.par
@@ -0,0 +1,9 @@
+# COINCOR parameter file
+
+input,s,a,,,,Input image root file name
+records,s,a,,,,Range of spectral records
+output,s,a,,,,Output file root name for new spectra
+start_rec,i,a,1,0,9999,Next starting spectral record
+ccmode,s,h,)_.ccmode,,,Correction mode (photo|iids)
+deadtime,r,h,)_.deadtime,,,Deadtime in seconds
+power,r,h,)_.power,,,IIDS power law coefficient
diff --git a/noao/onedspec/irsiids/coincor.x b/noao/onedspec/irsiids/coincor.x
new file mode 100644
index 00000000..df572c02
--- /dev/null
+++ b/noao/onedspec/irsiids/coincor.x
@@ -0,0 +1,123 @@
+# Coincidence correction options
+define	CC_PHOTO_MODE	1	# Photometer style correction
+define	CC_IIDS_MODE	2	# IIDS style
+define	CC_POWER_MODE	3	# Power law correction
+define	CC_USER_MODE	4	# User supplies a function
+
+
+# COINCOR -- Coincidence correction for detector deadtime
+
+procedure coincor (input, output, npts, expo, coflag, dt, power, mode)
+
+real	input[npts]
+real	output[npts]
+real	expo
+int	coflag
+real	dt
+real	power
+int	mode, npts
+
+begin
+	# Check that exposure time is legit
+	if (expo <= 0.0)
+	    return
+
+	# Select the method by which the correction is performed
+	switch (mode) {
+	    case CC_PHOTO_MODE:
+	    	# Photoelectric photometer
+		call ccphoto (input, output, npts, expo, coflag, dt)
+
+	    case CC_IIDS_MODE:
+		# IIDS style correction
+		if (coflag == -1) {
+		    call cciids (input, output, npts, expo, coflag, dt)
+		    if (power != 1.0)
+		        call ccpower (output, output, npts, expo, coflag,
+			    power)
+		} else if ((coflag == 0) && (power != 1.0))
+		    call ccpower (input, output, npts, expo, coflag, power)
+		else
+		    call amovr (input, output, npts)
+
+	    case CC_USER_MODE:
+		# Provided by the user
+		call ccuser (input, output, npts, expo, coflag, dt)
+	}
+end
+
+# CCPHOTO -- Photoelectric photometer coincidence correction
+
+procedure ccphoto (input, output, npts, expo, coflag, dt)
+
+real	input[npts], output[npts], expo, dt
+int	coflag
+int	npts
+
+int	i
+
+begin
+	do i = 1, npts
+	    output[i] = input[i] * exp (input[i] * dt / expo)
+	coflag = 2
+end
+
+# CCUSER -- User supplied correction scheme
+
+procedure ccuser (input, output, npts, expo, coflag, dt)
+
+real	input[npts], output[npts], expo, dt
+int	coflag
+int	npts
+
+begin
+	coflag = 3
+end
+
+# CCIIDS -- IIDS style correction scheme
+# From Instrumentation for Astronomy III (SPIE Vol 172) p.88 by Larry Goad
+#
+# Note that only the "Detect" mode of observation is supported.
+
+procedure cciids (input, output, npts, expo, coflag, dt)
+
+real	input[npts], output[npts], expo, dt
+int	npts, coflag
+
+int	i
+real	tsweep, value
+
+begin
+	# Allow tsweep to be the deadtime so that a different value
+	# may be entered for other instruments.
+	# For the IIDS, tsweep = 1.424e-3 sec
+	tsweep = dt
+
+	do i = 1, npts {
+	    value = 1 - input[i] / expo * tsweep
+	    if ((value < 0.) || (value > 1.))
+		output[i] = input[i]
+	    else
+	        output[i] = -expo * log (value)/ tsweep
+	}
+	coflag = 0
+end
+
+# CCPOWER -- Power law correction
+# Power law correction from Massey and De Veny, NOAO Newsletter #6.
+
+procedure ccpower (input, output, npts, expo, coflag, power)
+
+real	input[npts], output[npts], expo, power
+int	npts, coflag
+
+int	i
+
+begin
+	do i = 1, npts
+	    if (input[i] > 0.)
+	        output[i] = expo * (input[i] / expo) ** power
+	    else
+		output[i] = input[i]
+	coflag = 1
+end
diff --git a/noao/onedspec/irsiids/conversion.x b/noao/onedspec/irsiids/conversion.x
new file mode 100644
index 00000000..9d0b8c15
--- /dev/null
+++ b/noao/onedspec/irsiids/conversion.x
@@ -0,0 +1,213 @@
+define	MAX_CHARS	256
+
+
+# ASCII_TO_EBCDIC --  Vector procedure to convert ASCII characters to EBCDIC
+# characters using the lookup table atoe.
+
+procedure ascii_to_ebcdic (inbuffer, outbuffer, nchars)
+
+char	inbuffer[ARB]
+short	outbuffer[ARB], atoe[MAX_CHARS]
+int	l, nchars
+
+data (atoe[l], l = 1, 8)      / 0b,   1b,   2b,   3b, '7' , '-' , '.' , '/'  / 
+data (atoe[l], l = 9, 16)     /26b,   5b, '%' ,  13b,  14b,  15b,  16b,  17b /
+data (atoe[l], l = 17, 24)    /20b,  21b,  22b,  23b, '<' , '=' , '2' , '&'  / 
+data (atoe[l], l = 25, 32)    /30b,  31b, '?' , '\'',  34b,  35b,  36b,  37b / 
+data (atoe[l], l = 33, 40)    /'@' , 'O' , 177b, '{' , '[' , 'l' , 'P' , '}' / 
+data (atoe[l], l = 41, 48)    /'M' , ']' , '\\' , 'N' , 'k' , '`' , 'K' , 'a'/ 
+data (atoe[l], l = 49, 56)    /360b, 361b, 362b, 363b, 364b, 365b, 366b, 367b/ 
+data (atoe[l], l = 57, 64)    /370b, 371b, 'z' , '^' , 'L' , '~' , 'n' , 'o' / 
+data (atoe[l], l = 65, 72)    /'|' , 301b, 302b, 303b, 304b, 305b, 306b, 307b/ 
+data (atoe[l], l = 73, 80)    /310b, 311b, 321b, 322b, 323b, 324b, 325b, 326b/ 
+data (atoe[l], l = 81, 88)    /327b, 330b, 331b, 342b, 343b, 344b, 345b, 346b/ 
+data (atoe[l], l = 89, 96)    /347b, 350b, 351b, 'J' , 340b, 'Z' , '_' , 'm' / 
+data (atoe[l], l = 97, 104)   /'y' , 201b, 202b, 203b, 204b, 205b, 206b, 207b/ 
+data (atoe[l], l = 105, 112)  /210b, 211b, 221b, 222b, 223b, 224b, 225b, 226b/ 
+data (atoe[l], l = 113, 120)  /227b, 230b, 231b, 242b, 243b, 244b, 245b, 246b/ 
+data (atoe[l], l = 121, 128)  /247b, 250b, 251b, 300b, 'j' , 320b, 241b,   7b/ 
+data (atoe[l], l = 129, 136)  /' ' , '!' , '"' , '#' , '$' ,  25b,   6b,  27b/ 
+data (atoe[l], l = 137, 144)  /'(' , ')' , '*' , '+' , ',' ,  11b,  12b,  33b/ 
+data (atoe[l], l = 145, 152)  /'0' , '1' ,  32b, '3' , '4' , '5' , '6' ,  10b/ 
+data (atoe[l], l = 153, 160)  /'8' , '9' , ':' , ';' ,   4b,  24b, '>' , 341b/ 
+data (atoe[l], l = 161, 168)  /'A' , 'B' , 'C' , 'D' , 'E' , 'F' , 'G' , 'H' / 
+data (atoe[l], l = 169, 176)  /'I' , 'Q' , 'R' , 'S' , 'T' , 'U' , 'V' , 'W' / 
+data (atoe[l], l = 177, 184)  /'X' , 'Y' , 'b' , 'c' , 'd' , 'e' , 'f' , 'g' / 
+data (atoe[l], l = 185, 192)  /'h' , 'i' , 'p' , 'q' , 'r' , 's' , 't' , 'u' / 
+data (atoe[l], l = 193, 200)  /'v' , 'w' , 'x' , 200b, 212b, 213b, 214b, 215b/ 
+data (atoe[l], l = 201, 208)  /216b, 217b, 220b, 232b, 233b, 234b, 235b, 236b/ 
+data (atoe[l], l = 209, 216)  /237b, 240b, 252b, 253b, 254b, 255b, 256b, 257b/ 
+data (atoe[l], l = 217, 224)  /260b, 261b, 262b, 263b, 264b, 265b, 266b, 267b/ 
+data (atoe[l], l = 225, 232)  /270b, 271b, 272b, 273b, 274b, 275b, 276b, 277b/ 
+data (atoe[l], l = 233, 240)  /312b, 313b, 314b, 315b, 316b, 317b, 332b, 333b/ 
+data (atoe[l], l = 241, 248)  /334b, 335b, 336b, 337b, 352b, 353b, 354b, 355b/ 
+data (atoe[l], l = 249, 256)  /356b, 357b, 372b, 373b, 374b, 375b, 376b, 377b/
+
+begin
+	call alutcs (inbuffer, outbuffer, nchars, atoe)
+end
+
+# EBCDIC_TO_ASCII -- Vector procedure to convert EBCDIC characters to ASCII
+# characters.
+
+procedure ebcdic_to_ascii (inbuffer, outbuffer, nchars)
+
+char	outbuffer[ARB]
+short	inbuffer[ARB], etoa[MAX_CHARS]
+int	l, nchars
+
+data (etoa[l], l = 1, 8)      / 0b,   1b,   2b,   3b, 234b,  11b, 206b, 177b / 
+data (etoa[l], l = 9, 16)     /227b, 215b, 216b,  13b,  14b,  15b,  16b,  17b/ 
+data (etoa[l], l = 17, 24)    /20b,  21b,  22b,  23b, 235b, 205b,  10b, 207b / 
+data (etoa[l], l = 25, 32)    /30b,  31b, 222b, 217b,  34b,  35b,  36b,  37b / 
+data (etoa[l], l = 33, 40)    /200b, 201b, 202b, 203b, 204b,  12b,  27b,  33b/ 
+data (etoa[l], l = 41, 48)    /210b, 211b, 212b, 213b, 214b,   5b,   6b,   7b/ 
+data (etoa[l], l = 49, 56)    /220b, 221b,  26b, 223b, 224b, 225b, 226b,   4b/ 
+data (etoa[l], l = 57, 64)    /230b, 231b, 232b, 233b,  24b,  25b, 236b,  32b/ 
+data (etoa[l], l = 65, 72)    /' ' , 240b, 241b, 242b, 243b, 244b, 245b, 246b/ 
+data (etoa[l], l = 73, 80)    /247b, 250b, '[' , '.' , '<' , '(' , '+' , '!' / 
+data (etoa[l], l = 81, 88)    /'&' , 251b, 252b, 253b, 254b, 255b, 256b, 257b/ 
+data (etoa[l], l = 89, 96)    /260b, 261b, ']' , '$' , '*' , ')' , ';' , '^' / 
+data (etoa[l], l = 97, 104)   /'-' , '/' , 262b, 263b, 264b, 265b, 266b, 267b/ 
+data (etoa[l], l = 105, 112)  /270b, 271b, '|' , ',' , '%' , '_' , '>' , '?' / 
+data (etoa[l], l = 113, 120)  /272b, 273b, 274b, 275b, 276b, 277b, 300b, 301b/ 
+data (etoa[l], l = 121, 128)  /302b, '`' , ':' , '#' , '@' , '\'' , '=' , '"'/ 
+data (etoa[l], l = 129, 136)  /303b, 'a' , 'b' , 'c' , 'd' , 'e' , 'f' , 'g' / 
+data (etoa[l], l = 137, 144)  /'h' , 'i' , 304b, 305b, 306b, 307b, 310b, 311b/ 
+data (etoa[l], l = 145, 152)  /312b, 'j' , 'k' , 'l' , 'm' , 'n' , 'o' , 'p' / 
+data (etoa[l], l = 153, 160)  /'q' , 'r' , 313b, 314b, 315b, 316b, 317b, 320b/ 
+data (etoa[l], l = 161, 168)  /321b, '~' , 's' , 't' , 'u' , 'v' , 'w' , 'x' / 
+data (etoa[l], l = 169, 176)  /'y' , 'z' , 322b, 323b, 324b, 325b, 326b, 327b/ 
+data (etoa[l], l = 177, 184)  /330b, 331b, 332b, 333b, 334b, 335b, 336b, 337b/ 
+data (etoa[l], l = 185, 192)  /340b, 341b, 342b, 343b, 344b, 345b, 346b, 347b/ 
+data (etoa[l], l = 193, 200)  /'{' , 'A' , 'B' , 'C' , 'D' , 'E' , 'F' , 'G' / 
+data (etoa[l], l = 201, 208)  /'H' , 'I' , 350b, 351b, 352b, 353b, 354b, 355b/ 
+data (etoa[l], l = 209, 216)  /'}' , 'J' , 'K' , 'L' , 'M' , 'N' , 'O' , 'P' / 
+data (etoa[l], l = 217, 224)  /'Q' , 'R' , 356b, 357b, 360b, 361b, 362b, 363b/ 
+data (etoa[l], l = 225, 232)  /'\\', 237b, 'S' , 'T' , 'U' , 'V' , 'W' , 'X' / 
+data (etoa[l], l = 233, 240)  /'Y' , 'Z' , 364b, 365b, 366b, 367b, 370b, 371b/ 
+data (etoa[l], l = 241, 248)  /'0' , '1' , '2' , '3' , '4' , '5' , '6' , '7' / 
+data (etoa[l], l = 249, 256)  /'8' , '9' , 372b, 373b, 374b, 375b, 376b, 377b/
+
+begin
+	call alutsc (inbuffer, outbuffer, nchars, etoa)
+end
+
+# IBM_TO_ASCII -- Vector procedure for converting IBM characters to ASCII
+# characters.
+
+procedure ibm_to_ascii (inbuffer, outbuffer, nchars)
+
+char	outbuffer[ARB]
+short	inbuffer[ARB], ibmtoa[MAX_CHARS]
+int	l, nchars
+
+data (ibmtoa[l], l = 1, 8)     /0b,   1b,   2b,   3b, 234b,  11b, 206b, 177b  /
+data (ibmtoa[l], l = 9, 16)    /1227b, 215b, 216b,  13b,  14b,  15b,  16b, 17b/ 
+data (ibmtoa[l], l = 17, 24)   /20b,  21b,  22b,  23b, 235b, 205b,  10b, 207b / 
+data (ibmtoa[l], l = 25, 32)   /30b,  31b, 222b, 217b,  34b,  35b,  36b,  37b / 
+data (ibmtoa[l], l = 33, 40)   /200b, 201b, 202b, 203b, 204b,  12b,  27b,  33b/ 
+data (ibmtoa[l], l = 41, 48)   /210b, 211b, 212b, 213b, 214b,   5b,   6b,   7b/ 
+data (ibmtoa[l], l = 49, 56)   /220b, 221b,  26b, 223b, 224b, 225b, 226b,   4b/ 
+data (ibmtoa[l], l = 57, 64)   /230b, 231b, 232b, 233b,  24b,  25b, 236b,  32b/ 
+data (ibmtoa[l], l = 65, 72)   /' ' , 240b, 241b, 242b, 243b, 244b, 245b, 246b/ 
+data (ibmtoa[l], l = 73, 80)   /247b, 250b,   0b, '.' , '<' , '(' , '+' , '|' / 
+data (ibmtoa[l], l = 81, 88)   /'&' , 251b, 252b, 253b, 254b, 255b, 256b, 257b/ 
+data (ibmtoa[l], l = 89, 96)   /260b, 261b, '!' , '$' , '*' , ')' , ';' , '^' / 
+data (ibmtoa[l], l = 97, 104)  /'-' , '/' , 262b, 263b, 264b, 265b, 266b, 267b/ 
+data (ibmtoa[l], l = 105,112)  /270b, 271b,   0b, ',' , '%' , '_' , '>' , '?' / 
+data (ibmtoa[l], l = 113, 120) /272b, 273b, 274b, 275b, 276b, 277b, 300b, 301b/ 
+data (ibmtoa[l], l = 121, 128) /302b, '`' , ':' , '#' , '@' , '\'' , '=' , '"'/ 
+data (ibmtoa[l], l = 129, 136) /303b, 'a' , 'b' , 'c' , 'd' , 'e' , 'f' , 'g' / 
+data (ibmtoa[l], l = 137, 144) /'h' , 'i' , 304b, 305b, 306b, 307b, 310b, 311b/ 
+data (ibmtoa[l], l = 145, 152) /312b, 'j' , 'k' , 'l' , 'm' , 'n' , 'o' , 'p' / 
+data (ibmtoa[l], l = 153, 160) /'q' , 'r' , 313b, 314b, 315b, 316b, 317b, 320b/ 
+data (ibmtoa[l], l = 161, 168) /321b, '~' , 's' , 't' , 'u' , 'v' , 'w' , 'x' / 
+data (ibmtoa[l], l = 169, 176) /'y' , 'z' , 322b, 323b, 324b, 325b, 326b, 327b/ 
+data (ibmtoa[l], l = 177, 184) /330b, 331b, 332b, 333b, 334b, 335b, 336b, 337b/ 
+data (ibmtoa[l], l = 185, 192) /340b, 341b, 342b, 343b, 344b, 345b, 346b, 347b/ 
+data (ibmtoa[l], l = 193, 200) /'{' , 'A' , 'B' , 'C' , 'D' , 'E' , 'F' , 'G' / 
+data (ibmtoa[l], l = 201, 208) /'H' , 'I' , 350b, 351b, 352b, 353b, 354b, 355b/ 
+data (ibmtoa[l], l = 209, 216) /'}' , 'J' , 'K' , 'L' , 'M' , 'N' , 'O' , 'P' / 
+data (ibmtoa[l], l = 217, 224) /'Q' , 'R' , 356b, 357b, 360b, 361b, 362b, 363b/ 
+data (ibmtoa[l], l = 225, 232) /'\\', 237b, 'S' , 'T' , 'U' , 'V' , 'W' , 'X' / 
+data (ibmtoa[l], l = 233, 240) /'Y' , 'Z' , 364b, 365b, 366b, 367b, 370b, 371b/ 
+data (ibmtoa[l], l = 241, 248) /'0' , '1' , '2' , '3' , '4' , '5' , '6' , '7' / 
+data (ibmtoa[l], l = 249, 256) /'8' , '9' , 372b, 373b, 374b, 375b, 376b, 377b /
+
+begin
+	call alutsc (inbuffer, outbuffer, nchars, ibmtoa)
+end
+
+# Vector procedure to convert ASCII characters to ibm characters
+
+procedure ascii_to_ibm (inbuffer, outbuffer, nchars)
+
+char	inbuffer[ARB]
+short	outbuffer[ARB], atoibm[MAX_CHARS]
+int	l, nchars
+
+data (atoibm[l], l = 1, 8)     /0b,   1b,   2b,   3b, '7' , '-' , '.' , '/'   / 
+data (atoibm[l], l = 9, 16)    /26b,   5b, '%' ,  13b,  14b,  15b,  16b,  17b / 
+data (atoibm[l], l = 17, 24)   /20b,  21b,  22b,  23b, '<' , '=' , '2' , '&'  / 
+data (atoibm[l], l = 25, 32)   /30b,  31b, '?' , '\'',  34b,  35b,  36b,  37b / 
+data (atoibm[l], l = 33, 40)   /'@' , 'Z' , 177b, '{' , '[' , 'l' , 'P' , '}' / 
+data (atoibm[l], l = 41, 48)   /'M' , ']' , '\\', 'N' , 'k' , '`' , 'K' , 'a' / 
+data (atoibm[l], l = 49, 56)   /360b, 361b, 362b, 363b, 364b, 365b, 366b, 367b/ 
+data (atoibm[l], l = 57, 64)   /370b, 371b, 'z' , '^' , 'L' , '~' , 'n' , 'o' / 
+data (atoibm[l], l = 65, 72)   /'|' , 301b, 302b, 303b, 304b, 305b, 306b, 307b/ 
+data (atoibm[l], l = 73, 80)   /310b, 311b, 321b, 322b, 323b, 324b, 325b, 326b/ 
+data (atoibm[l], l = 81, 88)   /327b, 330b, 331b, 342b, 343b, 344b, 345b, 346b/ 
+data (atoibm[l], l = 89, 96)   /347b, 350b, 351b, 255b, 340b, 275b, '_' , 'm' / 
+data (atoibm[l], l = 97, 104)  /'y' , 201b, 202b, 203b, 204b, 205b, 206b, 207b/ 
+data (atoibm[l], l = 105, 112) /210b, 211b, 221b, 222b, 223b, 224b, 225b, 226b/ 
+data (atoibm[l], l = 113, 120) /227b, 230b, 231b, 242b, 243b, 244b, 245b, 246b/ 
+data (atoibm[l], l = 121, 128) /247b, 250b, 251b, 300b, 'O' , 320b, 241b,   7b/ 
+data (atoibm[l], l = 129, 136) /' ' , '!' , '"' , '#' , '$' ,  25b,   6b,  27b/ 
+data (atoibm[l], l = 137, 144) /'(' , ')' , '*' , '+' , ',' ,  11b,  12b,  33b/ 
+data (atoibm[l], l = 145, 152) /'0' , '1' ,  32b, '3' , '4' , '5' , '6' ,  10b/ 
+data (atoibm[l], l = 153, 160) /'8' , '9' , ':' , ';' ,   4b,  24b, '>' , 341b/ 
+data (atoibm[l], l = 161, 168) /'A' , 'B' , 'C' , 'D' , 'E' , 'F' , 'G' , 'H' / 
+data (atoibm[l], l = 169, 176) /'I' , 'Q' , 'R' , 'S' , 'T' , 'U' , 'V' , 'W' / 
+data (atoibm[l], l = 177, 184) /'X' , 'Y' , 'b' , 'c' , 'd' , 'e' , 'f' , 'g' / 
+data (atoibm[l], l = 185, 192) /'h' , 'i' , 'p' , 'q' , 'r' , 's' , 't' , 'u' / 
+data (atoibm[l], l = 193, 200) /'v' , 'w' , 'x' , 200b, 212b, 213b, 214b, 215b/ 
+data (atoibm[l], l = 201, 208) /216b, 217b, 220b, 232b, 233b, 234b, 235b, 236b/ 
+data (atoibm[l], l = 209, 216) /237b, 240b, 252b, 253b, 254b, 255b, 256b, 257b/ 
+data (atoibm[l], l = 217, 224) /260b, 261b, 262b, 263b, 264b, 265b, 266b, 267b/ 
+data (atoibm[l], l = 225, 232) /270b, 271b, 272b, 273b, 274b, 275b, 276b, 277b/ 
+data (atoibm[l], l = 233, 240) /312b, 313b, 314b, 315b, 316b, 317b, 332b, 333b/ 
+data (atoibm[l], l = 241, 248) /334b, 335b, 336b, 337b, 352b, 353b, 354b, 355b/ 
+data (atoibm[l], l = 249, 256) /356b, 357b, 372b, 373b, 374b, 375b, 376b, 377b/
+
+begin
+	call alutcs (inbuffer, outbuffer, nchars, atoibm)
+end
+
+# ALUTSC -- Vector operator to map one set of characters to another using a
+# lookup table.
+
+procedure alutsc (a, b, nchars, lut)
+
+char	b[nchars]
+int	nchars, i
+short	a[nchars], lut[ARB]
+
+begin
+	do i = 1, nchars, 1
+	    b[i] = lut[a[i] + 1]
+end
+
+# ALUTCS -- Vector operator to map one set of characters to another using
+# a lookup table.
+
+procedure alutcs (a, b, nchars, lut)
+
+char	a[nchars]
+int	nchars, i
+short	b[nchars], lut[ARB]
+
+begin
+	do i = nchars, 1, -1
+	    b[i] = lut[a[i] + 1]
+end
diff --git a/noao/onedspec/irsiids/doc/addsets.hlp b/noao/onedspec/irsiids/doc/addsets.hlp
new file mode 100644
index 00000000..6ce49122
--- /dev/null
+++ b/noao/onedspec/irsiids/doc/addsets.hlp
@@ -0,0 +1,66 @@
+.help addsets Feb85 noao.imred.iids/noao.imred.irs
+.ih
+NAME
+addsets - Add subsets of a string of spectra
+.ih
+USAGE
+addsets input records
+.ih
+PARAMETERS
+.ls input
+The root file name for the input spectra in the string.
+.le
+.ls records
+The range of spectra indicating the elements of the string.
+The names of the spectra will be formed by appending the range
+elements to the input root name.
+.le
+.ls output
+This is the root file name for the names of the spectra which will
+be created by the addset operation.
+.le
+.ls start_rec = 1
+The starting record number to be appended to the root name of the
+created spectra.
+.le
+.ls subset = 2
+The length of the substring of spectra which will be added together.
+For IIDS/IRS data which has been processed through BSWITCH, this
+parameter should be 2. This implies that spectra will be taken 
+2 at a time, added, and the sum written as a new spectrum.
+.le
+.ls weighting = yes
+If set to yes, an average of the substring of spectra is generated
+(if flux calibrated) weighted by the integration times of the
+individual spectra. If set to no, a simple average is generated.
+If not flux calibrated, this parameter has no effect - a simple
+sum is generated.
+.le
+.ih
+DESCRIPTION
+Every "subset" group of spectra will be accumulated and the sum will be
+written as a new spectrum. For example, if the input string contains
+100 spectra, and subset=2, then 50 new spectra will be created. Each
+new spectrum will be the sum of the consecutive pairs in the original string.
+
+If there are insufficient spectra to complete a subset accumulation,
+the sum is written out anyway and a warning printed. For example,
+if the input string contains 23 spectra, and subset=4, there will be
+6 new spectra created, but the last one will be based on only 3 spectra.
+
+Subset may be set to 1 to allow a copy operation although this is not
+a very efficient way to do so.
+.ih
+EXAMPLES
+The following three examples are those described above.
+
+.nf
+	cl> addsets nite1 2001-2100
+	cl> addsets nite1 2001-2023 subset=4
+	cl> addsets nite1 2001-2010 subset=1 output=nite2 \
+	>>> start_rec=2001
+.fi
+.ih
+SEE ALSO
+bswitch
+.endhelp
diff --git a/noao/onedspec/irsiids/doc/batchred.hlp b/noao/onedspec/irsiids/doc/batchred.hlp
new file mode 100644
index 00000000..9301f8b0
--- /dev/null
+++ b/noao/onedspec/irsiids/doc/batchred.hlp
@@ -0,0 +1,145 @@
+.help batchred Feb85 noao.imred.iids/noao.imred.irs
+.ih
+NAME
+batchred - Automated processing of IIDS/IRS spectra
+.ih
+USAGE
+batchred
+.ih
+PARAMETERS
+This script task has many parameters, but most are used as
+variables internal to the task and are not user parameters.
+There are 5 parameters having similar purposes: standard,
+sensfunc, bswitch, calibrate, and addsets. Each corresponds
+to the ONEDSPEC task of the same name and BATCHRED will generate
+the commands necessary to invoke those tasks if the associated
+parameter is set to yes (the default in all cases).
+
+.ls standard = yes
+.le
+.ls sensfunc = yes
+.le
+.ls bswitch = yes
+.le
+.ls calibrate = yes
+.le
+.ls addsets = yes
+.le
+.ls fnu = no
+This parameter is identical to the fnu parameter for CALIBRATE.
+.le
+.ls wave1 = 0.0
+This parameter is identical to the wave1 parameter for BSWITCH.
+.le
+.ls wave2 = 0.0
+This parameter is identical to the wave2 parameter for BSWITCH.
+.le
+.ls subset = 32767
+This parameter is identical to the subset parameter for BSWITCH.
+.le
+.ih
+DESCRIPTION
+Through a question and answer session, a series of commands to
+ONEDSPEC is generated which are then processed as a batch job
+to reduce "typical" spectra from the IIDS and IRS spectrographs.
+
+By setting the appropriate hidden parameters, the user may
+"turn off" command generation for any of the possible tasks.
+
+A script task is generated having the name "process.cl" which is
+submitted to the CL as the final command of BATCHRED.
+All terminal output which would normally appear during the course
+of running each of the individual tasks is redirected to a log file
+(default=ttylog).
+
+After the script has been generated, the user may suppress running
+the processing task. The script file remains on disk so that subsequent
+cases may be appended, such as when
+several independent runs of data are to be processed in one
+stream (e.g. several nights of data, each to be reduced separately).
+
+The questions which are asked are described below:
+
+"Root name for spectra file names:" This is the input root file name
+for all spectra which will be run through STANDARD and BSWITCH.
+
+"Root name for spectra to be created:" This is the output root file
+name which all newly created spectra will use. It is also the
+input file name for tasks CALIBRATE and ADDSETS since these tasks
+operate on spectra created by BSWITCH.
+
+"Starting record number for spectra to be created:" All created spectra
+will have a suffix number starting with this value and incremented
+by one for each new spectrum created.
+
+"File name to contain statistics information:" This file will contain
+informative output from SENSFUNC and BSWITCH. (default=stats)
+
+"File name to contain a log of terminal output:" All tasks talk back
+to let you know how things are proceding. The backtalk is saved
+in this file. (default=ttylog)
+
+"File name for output from STANDARD and input to SENSFUNC:" Just
+what it says. (default=std)
+
+"Record string to process:" The spectra are assumed to be representable
+by strings (try "help ranges" for details on the formats allowed).
+Both STANDARD and BSWITCH expect ranges of spectral record numbers
+which are appended to the root given in answer to the first question
+above. This question is asked repeatedly so that you can enter as
+many strings of spectra as you like and is ended by hitting return
+without entering a value. There is a short delay after entering
+each string of records while a check is made to verify that all
+your spectra actually exist.
+
+"Standard star name:" For each record string STANDARD expects
+the name of the standard star observed, but it must be given in
+a manner acceptable to STANDARD. (see STANDARD and LCALIB for
+more details).
+
+"Use weighted averages:" If answered yes, then SENSFUNC and BSWITCH
+will use their weighted averaging schemes.
+
+"Apply magnitude fudging:" If answered yes, then SENSFUNC will 
+use its "fudge" option. (see SENSFUNC)
+
+"Solve for grey additive extinction constant:" If answered yes, then
+SENSFUNC will solve for this value.
+
+"File name for sensitivity image file:" This will be the root name
+for the output sensitivity spectra from SENSFUNC.
+
+At anytime during the processing phase, you can inquire about the
+progress by listing the latest contents of the file "ttylog"
+either by "type ttylog" or by "tail ttylog". The latter command
+lists the last 12 lines of the file.
+
+Be sure to have all your record strings, standard star names,
+and options well planned and written down so that you can enter
+the answers correctly. The batch reductions are not overly
+tolerant of incorrect entries although some preliminary checks
+are performed during the entry process.
+
+.ih
+EXAMPLES
+
+The following invokes the batch reductions using all task options;
+
+	cl> batchred
+
+The following inhibits the STANDARD and SENSFUNC tasks which must have
+been run previously. This is equivalent to the IPPS "autoreduce":
+
+	cl> batchred standard- sensfunc-
+.ih
+BUGS
+If you make an error while entering the requested information, there
+is no way to effect repairs other than to (1) start all over, or (2) edit
+the generated script file "process.cl" using the system editor.
+
+If a task encounters an irrecoverable error, the background job
+hangs until you kill it using "kill N" where N is the job number.
+.ih
+SEE ALSO
+mkscript, standard, sensfunc, bswitch, calibrate, addsets
+.endhelp
diff --git a/noao/onedspec/irsiids/doc/bswitch.hlp b/noao/onedspec/irsiids/doc/bswitch.hlp
new file mode 100644
index 00000000..a50647b4
--- /dev/null
+++ b/noao/onedspec/irsiids/doc/bswitch.hlp
@@ -0,0 +1,228 @@
+.help bswitch Sep87 noao.imred.iids/noao.imred.irs
+.ih
+NAME
+bswitch - generate sky-subtracted accumulated spectra
+.ih
+USAGE
+bswitch input records
+.ih
+PARAMETERS
+.ls input
+The root name for the input spectra to be beam-switched.
+.le
+.ls records
+The range of spectra to be included in the beam-switch operation.
+Each range item will be appended to the root name to form an image
+name. For example, if "input" is "nite1" and records is "1011-1018",
+then spectra nite1.1011, nite.1012 ... nite1.1018 will be included.
+.le
+.ls output
+New spectra are created by the beam-switch operation. This parameter
+specifies the root name to be used for the created spectra.
+.le
+.ls start_rec = 1
+Each new spectrum created has "output" as its root name and a trailing
+number appended. The number begins with start_rec and is incremented
+for each new spectrum. For example, if "output" is given as "nite1b"
+and start_rec is given as 1001, then new spectra will be created as
+nite1b.1001, nite1b.1002 ...
+.le
+.ls stats = "stats"
+A file by this name will have statistical data appended to it, or created
+if necessary. If a null file name is given (""), no statistical output
+is given. For each aperture, a listing of countrates for each
+observation is given relative to the observation with the highest rate.
+.le
+.ls ids_mode = yes
+If the data are taken under the usual IIDS "beam-switch" mode, this
+parameter should be set to yes so that accumulations will be performed
+in pairs. But if the data are taken where there is no sky observation
+or different numbers of sky observations, ids_mode should be set to no.
+If weighting is in effect, ids_mode=yes implies weighting of the
+object-sky sum; if ids_mode=no, then weighting is applied to the
+object and sky independently because then there is no guarantee that
+an object and sky observation are related.
+.le
+.ls extinct = yes
+If set to yes, a correction for atmospheric extinction is applied.
+The image header must have either a valid entry for AIRMASS or
+for hour angle (or right ascension and sidereal time) and declination.
+.le
+.ls weighting = no
+If set to yes, the entire spectrum or a specified region will be used
+to obtain a countrate indicative of the statistical weight to be
+applied to the spectrum during the accumulations.
+.le
+.ls subset = 32767
+A subset value larger than the number of independent spectra to be
+added indicates that the operation is to produce a single spectrum
+for each aperture regardless of how many input spectra are entered.
+If subset is a smaller number, say 4, then the accumulations
+are written out after every 4 spectra and then re-initialized to zero
+for the next 4.
+.le
+.ls wave1 = 0.0
+If weighting=yes, this parameter indicates the starting point in the
+spectrum for the countrate to be assessed. For emission-line objects,
+this is particularly useful because the regime of information is then
+confined to a narrow spectral region rather than the entire spectrum.
+Defaults to the beginning of the spectrum.
+.le
+.ls wave2 = 0.0
+This provides the ending wavelength for the countrate determination.
+Defaults to the endpoint of the spectrum.
+.le
+.ls observatory = "observatory"
+Observatory at which the spectra were obtained if
+not specified in the image header by the keyword OBSERVAT.  The
+observatory may be one of the observatories in the observatory
+database, "observatory" to select the observatory defined by the
+environment variable "observatory" or the task \fBobservatory\fR, or
+"obspars" to select the current parameters set in the \fBobservatory\fR
+task.  See help for \fBobservatory\fR for additional information.
+.le
+.ls extinction = ")_.extinction"
+The the name of the file containing extinction values.
+Required if extinct=yes.
+.le
+.ih
+DESCRIPTION
+Data from multiaperture spectrographs are summed according to
+aperture number and sky subtracted if sky observations are available.
+Data for up to 50 apertures may be simultaneously accumulated.
+The accumulated spectra are written to new images. 
+
+The exposure times for each observation may be different. All
+internal computations are performed in terms of count rates,
+and converted back to counts (for statistical analysis) prior to writing
+the new image. Therefore, the time on the sky and object may
+be different as well. When these extensions to the normal
+mode are required, the flag ids_mode must be set to no.
+Then object and sky accumulations are performed totally
+independently and a difference is derived at the conclusion
+of the operation.
+
+If ids_mode is set to yes, then the usual IIDS/IRS "beam-switch"
+observing mode is assumed. This implies that an equal number of
+sky and object spectra are obtained through each aperture
+after 2N spectra have been accumulated, where N is the number
+of instrument apertures (2 for the IIDS/IRS). It is also assumed
+that the object and sky exposure times are equal for each aperture.
+Note that the "nebular" mode (where all instrument apertures
+point at an extended object simultaneously, and then all apertures
+point at sky simultaneously) is an acceptable form for
+beam-switched data in ids_mode.
+
+The accumulations are optionally weighted by the countrate
+over a region of the spectrum to improve the statistics during
+variable conditions. The user may specify the region of spectrum
+by wavelength. In ids_mode, the statistics are obtained from
+object-sky differences; otherwise, the statistics are performed
+on object+sky and sky spectra separately.
+
+The spectra may be extinction corrected if this has not already
+been performed.
+In order to perform either the extinction correction or the
+weighting process, the spectra must have been placed on a linear
+wavelength scale (or linear in the base 10 logarithm).
+
+Strings of spectra are  accumulated to produce a single
+summed spectrum for each observing aperture. But in some cases
+it is desirable to produce summed spectra from subsets of the
+entire string to evaluate the presence of variations either due
+to observing conditions or due to the physical nature of the
+object. A subset parameter may be set to the frequency at which
+spectra are to be summed.
+
+In order that the processing occur with minimal user interaction,
+elements from the extended image header are used to direct the
+flow of operation and to obtain key observing parameters.
+The required parameters are: object/sky flag (OFLAG=1/0), exposure
+time in seconds (ITM), beam (that is, aperture) number (BEAM-NUM), airmass (AIRMASS)
+or alternatively hour angle (HA) and declination (DEC), or
+right ascension (RA), sidereal time (ST), declination (DEC), and the
+observatory (OBSERVAT),
+starting wavelength (W0), and wavelength increment per channel (WPC),
+where the names in parenthesis are the expected keywords in the
+header.  If the observatory is not specified in the image the
+observatory parameter is used.  See \fBobservatory\fR for further
+details on the observatory database.
+
+The following header flags are used as well: DC_FLAG
+for dispersion corrected data (must=0), BS_FLAG for beam-switching
+(must not be 1 which indicates the operation was already done),
+EX_FLAG for extinction correction (if = 0 extinction is assumed already
+done).  
+
+The headers may be listed with the IMHEADER task, setting
+the parameter "long" = yes. The values for the parameters follow 
+the rules used for IIDS and IRS data.
+
+After the beam-switch operation, the newly created spectra will
+have header elements taken from the last object spectrum.
+A few parameters will be updated to reflect the operation
+(e.g. integration time, processing flags).
+
+.ih
+EXAMPLES
+The following example will accumulate a series of 16 spectra obtained
+in the normal beam-switched mode and create two new extinction corrected
+spectra having names nite1bs.1 and nite1bs.2:
+
+	cl> bswitch nite1 1011-1026 nite1bs 1
+
+The following example performs the same functions but accumulates the data
+to produce 8 new spectra representing the individual object-sky pairs:
+
+	cl> bswitch nite1 1011-1026 nite1bs 1 subset=4
+
+The following example produces an extinction corrected spectrum for every
+input spectrum. Note that ids_mode is set to off to generate separate object and
+sky sums, and subset is set to 2 so that every pair of spectra (one object and
+one sky) are written out as two new spectra:
+
+	cl> bswitch nite1 1011-1026 nite1bs 1 subset=2 ids_mode-
+
+The next example produces a pair of spectra for each of 3 independent
+objects observed, provided that each was observed for the same number
+of observations (16 in this case).
+
+.nf
+	cl> bswitch nite1 1011-1026,1051-1066,1081-1096 nite1bs 1 \
+	>>> subset=16
+.fi
+
+The next example shows how to use the weighting parameters where
+the indicative flux is derived from the region around the emission-line
+of 5007A.
+
+.nf
+	cl> bswitch nite1 1011-1026 nite1bs 1 weighting- \
+	>>> wave1=4990, wave2=5020
+.fi
+.ih
+TIME REQUIREMENTS
+The principle time expenditure goes toward extinction correcting the
+data. For IIDS type spectra (length=1024 pixels), approximately 30 cpu
+seconds are required to beam-switch a series of 16 spectra.
+.ih
+BUGS
+The number of apertures is restricted to 50 and must be labeled
+between 0 and 49 in the image header (the IIDS uses 0 and 1).
+
+Until an image header editor is available, BSWITCH 
+can be applied only to data with properly prepared headers
+such as IIDS/IRS data read by RIDSMTN, RIDSFILE and some data via RFITS.
+
+When used to perform the function of extinction correction only (the
+third example above), the statistics file fails to note the output
+image name for the sky spectrum.
+
+The data must be on a linear wavelength scale.
+The starting wavelength, W0, and a wavelength
+per channel, WPC, are required header information, and the DC_FLAG
+must be set to 0.
+.ih
+SEE ALSO
+observatory, sensfunc, imheader, lcalib, ridsmtn, ridsfile, rfits
+.endhelp
diff --git a/noao/onedspec/irsiids/doc/coefs.hlp b/noao/onedspec/irsiids/doc/coefs.hlp
new file mode 100644
index 00000000..777933bc
--- /dev/null
+++ b/noao/onedspec/irsiids/doc/coefs.hlp
@@ -0,0 +1,57 @@
+.help coefs May85 noao.imred.iids/noao.imred.irs
+.ih
+NAME
+coefs -- Extract dispersion coefs from mtn HeNeAr headers
+.ih
+USAGE
+coefs input records database
+.ih
+PARAMETERS
+.ls input
+The input image root name for the spectral images containing the
+dispersion coefficients.
+.le
+.ls records
+The range of records for which the root name applies.
+.le
+.ls database
+The database file name which will contain the coefficients.
+.le
+.ih
+DESCRIPTION
+The spectra specified by the combination of the root name
+and the records are scanned for the presence of dispersion
+coefficients. If present, the coefficients and necessary
+information are written to the file indicated by the database
+parameter. This file an then be used by the linearization
+program DISPCOR to correct any spectra for which the
+database is appropriate.
+
+Each invocation of COEFS appends to the database file, or
+creates a new file if necessary.
+
+The following assumptions are made concerning the coefficients,
+which are always correct for IIDS and IRS mountain reduced
+data at Kitt Peak.
+.ls  5 (1)
+The coefficients represent Legendre polynomials.
+.le
+.ls (2)
+The coefficients apply to pixels 1 through 1024 in the original data.
+.le
+.ih
+EXAMPLES
+The following example reads the coefficients from the headers
+for nite1 arc spectra taken near the beginning and end of the
+night and creates a database file called nite1.db:
+
+	cl> coefs nite1 3-4,201-202 nite1.db
+
+.ih
+TIME REQUIREMENTS
+Approximately 1 second per spectrum is required. This is primarily
+overhead due to file access.
+.ih
+SEE ALSO
+dispcor, identify
+.endhelp
diff --git a/noao/onedspec/irsiids/doc/coincor.hlp b/noao/onedspec/irsiids/doc/coincor.hlp
new file mode 100644
index 00000000..74e002f3
--- /dev/null
+++ b/noao/onedspec/irsiids/doc/coincor.hlp
@@ -0,0 +1,101 @@
+.help coincor Feb87 noao.imred.iids/noao.imred.irs
+.ih
+NAME
+coincor -- Correct detector count rates
+.ih
+USAGE
+coincor input records
+.ih
+PARAMETERS
+.ls input
+The root file name of the input spectra.
+.le
+.ls records
+The range of spectra.
+The names of the spectra will be formed by appending the range
+elements to the input root name.
+.le
+.ls output
+This is the root file name for the corrected spectra.  If no root name
+is specified (specified with the null string "") then the operation
+is done in place.
+.le
+.ls start_rec = 1
+The starting record number to be appended to the root name of the
+created spectra.
+.le
+.ls ccmode = )_.ccmode
+The mode used to model the detector count rate corrections.
+In the following C(obs) is the observed count rate and C(cor) is the
+corrected count rate.
+.ls "photo"
+Photoelectric photometer with discriminator mode.  The count rate
+correction is
+
+    C(cor) = C(obs) * exp (C(obs) * deadtime)
+    
+where the parameter \fIdeadtime\fR is the representative deadtime in seconds.
+.le
+.ls "iids"
+IIDS correction given by
+
+    C(cor) = (-ln(1-C(obs)*deadtime)/deadtime)**power
+
+where \fBdeadtime\fR is a parameter related to the sweep time used to
+correct for coincidence losses and \fBpower\fR is a power law coefficient.
+.le
+.le
+.ls deadtime = )_.deadtime
+For the "photo" mode this parameter is the period, in seconds, during
+which no counts can be registered by the detector.  Note that this is
+based on a per pixel basis.  So if the discriminator dead period is of
+order 50 nanoseconds and 2000 pixels are observed per readout, the
+effective deadtime is about 10E-4 seconds.  For the "iids" mode this
+parameter defines the sweep time correction and has a value of 1.424E-3
+seconds.
+.le
+.ls power = )_.power
+The IIDS power law coefficient.  The standard value is 0.975.
+.le
+.ih
+DESCRIPTION
+The input spectra are corrected for detector count rate errors.  If no
+output root name is given then the operation is done in place.  The type
+of correction is specified by the parameter \fIccmode\fR.  The available
+modes are for a general photomultiplier with discriminator coincidence
+correction, and the NOAO IIDS.  The parameters for these modes are
+\fIdeadtime\fR and \fIpower\fR.  The exposure time, in seconds, is a
+required image header parameter (keyword = EXPOSURE).
+
+The default mode is for the NOAO IIDS.  The IIDS correction includes a
+power law correction for a nonlinear effect in the IIDS image tube chain
+which is not included by the mountain reduction software at the telescope.
+If the spectra have been coincidence corrected at the telescope
+then only the nonlinear power law correction is applied.
+
+The coincidence correction flag may take the values -1 for no correction,
+0 for the IIDS correction with \fIpower\fR = 1 (the correction
+applied by the mountain reduction software), 1 for the full IIDS
+correction, and 2 for the photomuliplier mode correction.
+.ih
+EXAMPLES
+The following example corrects a series of IIDS spectra:
+
+	cl> coincor nite1 1-250 output=nite1cc start_rec=1
+
+The following example corrects a series of spectra from the
+Lick ITS:
+
+.nf
+	cl> coincor its 1-250 output=itscc start=1 ccmode=photo \
+	>>> deadtime=2.4E-4 power=1
+.fi
+.ih
+TIME REQUIREMENTS
+\fBCoincor\fR requires approximately 1 second per spectrum of length 1024.
+.ih
+SEE ALSO
+.nf
+The \fBimred.iids\fR package is designed for reducing NOAO IIDS spectra.
+.fi
+.endhelp
diff --git a/noao/onedspec/irsiids/doc/extinct.hlp b/noao/onedspec/irsiids/doc/extinct.hlp
new file mode 100644
index 00000000..66aca3d6
--- /dev/null
+++ b/noao/onedspec/irsiids/doc/extinct.hlp
@@ -0,0 +1,49 @@
+.help extinct Apr85 noao.onedspec
+.ih
+NAME
+extinct -- Correct spectra for atmospheric extinction
+.ih
+USAGE
+extinct root records output
+.ih
+PARAMETERS
+.ls root
+The root name for the input spectra to be corrected.
+.le
+.ls records
+The range of spectra to be included in the extinction operation.
+.le
+.ls output
+The root name for the output corrected spectra
+.le
+.ls start_rec
+The starting record number for the output corrected spectra.
+.le
+.ls nr_aps = 2
+The number of instrument apertures for this data set.
+.le
+.ih
+DESCRIPTION
+The input spectra are corrected for atmospheric extinction. 
+EXTINCT redirects the spectra through the task BSWITCH so all
+procedures are identical to those described for that task.
+
+Because BSWITCH attempts to perform a beam-switch operation
+unless the subset parameter is equal to the number of
+instrument apertures (in which case beam-switching degenerates
+to a copy operation), the hidden parameter nr_aps should be set
+appropriately to the instrument. For IIDS and IRS data, this
+is 2.
+.ih
+EXAMPLES
+
+	cl> extinct  nite1 1001-1032 nite1ex
+.ih
+BUGS
+The input string of spectra must be ordered so that only
+one spectrum from each aperture is present among substrings
+of length nr_aps.
+.ih
+SEE ALSO
+bswitch
+.endhelp
diff --git a/noao/onedspec/irsiids/doc/flatdiv.hlp b/noao/onedspec/irsiids/doc/flatdiv.hlp
new file mode 100644
index 00000000..e6e8c22e
--- /dev/null
+++ b/noao/onedspec/irsiids/doc/flatdiv.hlp
@@ -0,0 +1,94 @@
+.help flatdiv Dec86 noao.imred.iids/noao.imred.irs
+.ih
+NAME
+flatdiv -- Divide spectra by flat field spectra
+.ih
+USAGE
+flatdiv input records
+.ih
+PARAMETERS
+.ls input
+The root file name for the input records to be divided.
+.le
+.ls records
+The range of spectra to be included in the divide operation.
+Each range item will be appended to the root name to form an
+image file name.
+.le
+.ls output
+New spectra are created by the flatdiv operation. This parameter
+specifies the root name to be used for the created spectra.
+.le
+.ls start_rec
+Each new spectrum created as "output" has its root name and a 
+trailing number appended starting with "start_rec". Subsequent
+output images will have an incremented trailing number.
+Note that even if an output image is not created because the input
+image has already been flattened or the input image is not found the
+output record number is still incremented.
+.le
+.ls flat_file
+The root name for the sensitivity spectra as produced by FLATFIT.
+Normally with multi-aperture instruments, FLATFIT will produce a
+spectrum appropriate to each aperture and the file name will have
+"flat_file" as the file name root and the aperture number appended.
+.le
+.ls coincor = )_.coincor
+If set to yes, coincidence correction is applied to the data during
+the division, and the following three parameters are required.
+For more about this correction see \fBcoincor\fR.
+.ls ccmode = )_.ccmode
+The mode by which the coincidence correction is to be performed.
+This may be "iids" or "photo".
+.le
+.ls deadtime = )_.deadtime
+The detector deadtime in seconds.
+.le
+.ls power = )_.power
+Power law IIDS non-linear correction exponent.
+.le
+.le
+.ih
+DESCRIPTION
+The input spectra are divided by the flat fields which are
+represented by spectra produced by FLATFIT.
+
+To avoid possible division by zero, any zeroes in the flat field
+spectra generated by FLATFIT are replaced by 1.0.
+
+The input spectra may optionally be corrected for coincidence losses.
+
+If the input and output spectra (after appending the record numbers) are
+the same then the division is performed in-place; i.e. the flattened spectra
+replace the original input spectra.
+Note that even if an output image is not created because the input
+image has already been flattened or the input image is not found the
+output record number is still incremented.  This is to insure that if
+in-place division is desired that the input and output names remain
+matched.
+.ih
+EXAMPLES
+The following example divides a series of spectra to produce 20 new
+spectra having names nite1.1221 ... nite1.1240.
+
+	cl> flatdiv nite1 1201-1220 nite1 1221
+
+The same spectra as above are simultaneously corrected for
+coincidence losses.
+
+	cl> flatdiv nite1 1201-1220 nite1 1221 coincor=yes
+
+The flattened spectra replace the unflattened spectra.
+
+	cl> flatdiv nite1 1201-1220 nite1 1201
+
+Note that the input record numbers must be contiguous and the starting
+output record number must be the same as the first input record number.
+.ih
+TIME REQUIREMENTS
+Approximately 1 second is required to correct a spectrum of length
+1024 points.
+.ih
+SEE ALSO
+coincor, flatfit
+.endhelp
diff --git a/noao/onedspec/irsiids/doc/flatfit.hlp b/noao/onedspec/irsiids/doc/flatfit.hlp
new file mode 100644
index 00000000..af84cb3c
--- /dev/null
+++ b/noao/onedspec/irsiids/doc/flatfit.hlp
@@ -0,0 +1,188 @@
+.help flatfit Dec86 noao.imred.iids/noao.imred.irs
+.ih
+NAME
+flatfit -- Sum and normalize flat field spectra
+.ih
+USAGE
+flatfit root records
+.ih
+PARAMETERS
+.ls root
+The root file name for the input names of the flat field
+spectra to be accumulated and fit for normalization.
+.le
+.ls records
+The range of spectra indicating the elements of the string.
+The names of the spectra will be formed by appending the range
+elements to the input root name.
+.le
+.ls output
+This is the root file name for the names of the spectra which will
+be created during normalization. The aperture number for the observation
+will be appended to the root in form "root.nnnn" where nnnn is the aperture
+number with leading 0's.
+.le
+.ls function = "chebyshev"
+The accumulated spectra are fit by this function type - either
+chebyshev or legendre polynomials, or spline3 or spline1 interpolators.
+.le
+.ls order = 4
+The order of the fit using the above function. This should generally be
+a low order fit to avoid introduction of high spatial frequency wiggles.
+.le
+.ls niter = 1
+The number of iterations to reject discrepant pixels upon initial
+startup of the solution.
+.le
+.ls lower = 2.0
+The number of sigmas for which data values less than this cutoff are
+rejected.
+.le
+.ls upper = 2.0
+The number of sigmas for which data values greater than this cutoff are
+rejected.
+.le
+.ls ngrow = 0
+The number of pixels on either side of a rejected pixel to also be rejected.
+.le
+.ls div_min = 1.0
+During the normalization process, a division by zero will produce
+this value as a result.
+.le
+.ls interact = yes
+If set to yes, graphical interaction with the normalization process
+is provided for at least the first aperture for which sums are available.
+If set to no, no interaction is provided.
+.le
+.ls all_interact = no
+If set to yes, then interaction will be provided for all apertures
+for which sums have been accumulated. If set to no then the parameter interact
+will determine if the first aperture data is to be interactive.
+.le
+.ls coincor = )_.coincor
+If set to yes, coincidence correction is applied to the data during
+the summation process, and the following three parameters are required.
+See \fBcoincor\fR for more about this correction.
+.ls ccmode = )_.ccmode
+The mode by which the coincidence correction is to be performed.
+This may be "iids" or "photo".
+.le
+.ls deadtime = )_.deadtime
+The detector deadtime in seconds.
+.le
+.ls power = )_.power
+Power law IIDS non-linear correction exponent.
+.le
+.le
+.ls cursor = ""
+Graphics cursor input.  When null the standard cursor is used otherwise
+the specified file is used.
+.le
+.ih
+DESCRIPTION
+The specified spectra are added by aperture number to produce
+summations which are then fit by a specified fitting function.
+The fitting function is then divided into the sum to produce a
+normalized (to 1.0) sum in which the low frequency spatial
+response has been removed.
+
+The resultant normalized images may then be divided into all other
+data to remove the pixel-to-pixel variations without introducing
+any color terms. The spectra may be used directly if they happen
+to be object spectra in which the low frequency response is to be
+removed.
+
+During the accumulation process the spectra may be corrected for
+coincidence losses if the detector is subject to the phenomenon.
+
+After accumulating all input spectra, the pixels in each sum are
+fit according to
+the specified function. If the interactive switches are set, then
+graphical interaction is made available. If only the interact parameter
+is set to yes, then only the data from the first aperture will
+be available for interaction. Data from subsequent apertures will
+be fit using the same parameters and number of iterations as the
+first. If the all_interact parameter is also
+set, then data from each aperture will be presented for interaction.
+
+At each step in the fit, pixels which are discrepant by more than
+"upper" sigmas above the fit, or "lower" sigmas below the fit, are
+rejected. The rejection process may be applied many times (iterations)
+to continue rejecting pixels. If the upper and lower sigmas are
+not equal, the resulting fit will be biased slightly above the mean
+(for lower < upper) or below the mean (upper < lower). This is useful
+when the spectrum being fit is that of a star having either absorption
+or emission lines.
+ 
+A display is presented of the sum and the fit through the data.
+A status line is printed containing the fit type, the order of
+the fit, the rms residual from the fit, and the number of data
+points in the fit after one iteration of rejection.
+
+The following cursor keystrokes are then active:
+.ls ?
+Clear the screen and display the active keystrokes
+.le
+.ls /
+Indicate active keystrokes on the status line
+.le
+.ls e
+Change plot mode to an error plot. This display is defined
+as the deviation from the fit divided by the data values [ (data - fit)/ data]
+at each pixel
+.le
+.ls f
+Change plot mode back to the fit through the data display
+.le
+.ls o
+Change the order of the fit.
+.le
+.ls l
+Change the lower rejection criterion (in units of sigma).
+.le
+.ls u
+Change the upper rejection criterion.
+.le
+.ls s
+Change both rejection criteria to the same value.
+.le
+.ls r
+Reinstate rejected pixels.
+.le
+.ls i
+Iterate one more time.
+.le
+.ls n
+Iterate several more times - the user is prompted for the count.
+.le
+.ls q
+Quit and accept the solution
+.le
+.ls <CR>
+RETURN is the same as 'q' but a confirmation request to exit must be
+answered as yes.
+.le
+
+All keystrokes but ?,/,e,f, and q force another iteration which will
+reject additional pixels. To fully inhibit pixel rejection, the sigmas
+should be set to a large value (e.g. 100).
+.ih
+EXAMPLES
+The following example will accumulate 8 spectra and fit the first
+aperture data interactively but not the second, and apply coincidence
+corrections to the sums. The upper and lower rejection criteria
+have been altered to bias the seventh order fit to a higher level.
+
+	cl> flatfit nite1 1-4,201-204 coin+ low=1.4 up=3 order=7
+.ih
+BUGS
+For some reason, the error plot is supposed to have a zero level line
+drawn, but none appears.
+
+As in most of the IRAF software, the order of a fit refers to the number
+of terms in the fit, so that a fit of order 1 implies a constant and order
+2 implies a linear fit.
+.ih
+SEE ALSO
+coincor, flatdiv
+.endhelp
diff --git a/noao/onedspec/irsiids/doc/powercor.hlp b/noao/onedspec/irsiids/doc/powercor.hlp
new file mode 100644
index 00000000..e1f9c70e
--- /dev/null
+++ b/noao/onedspec/irsiids/doc/powercor.hlp
@@ -0,0 +1,62 @@
+.help powercor Oct86 noao.imred.iids/noao.imred.irs
+.ih
+NAME
+powercor -- Apply power law correction to mountain reduced spectra
+.ih
+USAGE
+powercor input records
+.ih
+PARAMETERS
+.ls input
+The root file name of the input spectra.
+.le
+.ls records
+The range of spectra.
+The names of the spectra will be formed by appending the range
+elements to the input root name.
+.le
+.ls output
+This is the root file name for the corrected spectra.
+.le
+.ls start_rec = 1
+The starting record number to be appended to the root name of the
+created spectra.
+.le
+.ls power = )iids.power
+The power law coefficient.
+.le
+.ih
+DESCRIPTION
+A power law correction to the IIDS count rates is applied to the input
+spectra.  The mountain reduction software applies a coincidence correction
+to the observed IIDS count rates but does not correct for a nonlinear effect
+in the image tube chain.  This second correction takes the form of a
+power law
+
+	C(out) = C(in) ** power
+
+where C(in) is the input, coincidence corrected, count rate and C(out)
+is the corrected count rate.  The power is a parameter of the task
+which defaults to the \fBiids\fR package parameter set to the appropriate
+value for the IIDS.  The exposure time, in seconds, is a required
+image header parameter (keyword = EXPOSURE) used to convert the
+total counts to count rates.
+
+Note that if the original raw spectra are being reduced then the either
+\fBcoincor\fR or \fBpowercor\fR may be used to apply both the coincidence
+correction and the power law correction at the same time.  In other words,
+the tasks apply the coincidence correction if the coincidence flag (CO-FLAG) is
+-1 (uncorrected) and the power law correction alone if the flag is zero
+(coincidence corrected only).  The flag is 1 when both the coincidence and
+nonlinear correction have been applied.
+
+This task is a script calling \fBcoincor\fR with \fIccmode\fR = "iids".
+.ih
+EXAMPLES
+The following example corrects a series of IIDS spectra:
+
+	cl> powercor nite1 1-250 output=nite1cc start_rec=1
+.ih
+SEE ALSO
+coincor
+.endhelp
diff --git a/noao/onedspec/irsiids/doc/process.hlp b/noao/onedspec/irsiids/doc/process.hlp
new file mode 100644
index 00000000..5cedcde3
--- /dev/null
+++ b/noao/onedspec/irsiids/doc/process.hlp
@@ -0,0 +1,20 @@
+.help process Oct85 noao.imred.iids/noao.imred.irs
+.ih
+NAME
+process	-- A task generated by BATCHRED
+.ih
+USAGE
+process
+.ih
+DESCRIPTION
+The task \fBbatchred\fR creates a script called process.cl for batch
+reductions.  \fBBatchred\fR also has an option to automatically run
+this script.
+.ih
+EXAMPLES
+The task \fBbatchred\fR is run to setup a set of beam switching operations.
+It creates the script \fBprocess.cl\fR which the user runs as a background
+process as follows:
+
+	cl> process&
+.endhelp
diff --git a/noao/onedspec/irsiids/doc/slist1d.hlp b/noao/onedspec/irsiids/doc/slist1d.hlp
new file mode 100644
index 00000000..6c7d2702
--- /dev/null
+++ b/noao/onedspec/irsiids/doc/slist1d.hlp
@@ -0,0 +1,59 @@
+.help slist1d Jan92 noao.imred.irs/iids
+.ih
+NAME
+slist1d -- List spectral header information
+.ih
+USAGE
+slist1d input records
+.ih
+PARAMETERS
+.ls input
+The image root name for the spectra to be listed.
+.le
+.ls records
+The record string for the spectra to be listed. The records will be appended
+to the root name to form image names of the type "root.xxxx".
+.le
+.ls long_header = no
+If set to yes, then a complete listing of the header elements
+is given. If set to no, then a single line per spectrum is given which lists
+in the following order: the image name, object or sky spectrum, exposure
+time, spectrum length, and image title.
+.le
+.ih
+DESCRIPTION
+Each spectrum in the list implied by the root name and the record string
+is opened and the header is read. The pixel file is not accessed in order
+to save time. The header listing is directed to STDOUT and may be
+redirected for printing.
+
+A warning message is issued if
+a requested image is not found, but otherwise proceeds.
+.ih
+EXAMPLES
+The following example lists 8 spectral headers in long form on the printer:
+
+.nf
+	cl> slist1d nite1 1001-1008 | lprint
+.fi
+
+The next example lists the same spectral headers but in short form
+on the terminal
+
+.nf
+	cl> slist1d nite1 1001-1008 long-
+.fi
+.ih
+REVISIONS
+.ls SLIST1D V2.10
+This task is the same as V2.9 \fBslist\fR and applies only to the older
+IRS/IIDS record extension spectra.  In V2.10 \fBslist\fR
+has been revised for multiaperture spectra.
+.le
+.ih
+BUGS
+SLIST1D does not inform the user if the pixel file can or cannot be read.
+.ih
+SEE ALSO
+slist, imheader
+.endhelp
diff --git a/noao/onedspec/irsiids/doc/subsets.hlp b/noao/onedspec/irsiids/doc/subsets.hlp
new file mode 100644
index 00000000..a9f0ae68
--- /dev/null
+++ b/noao/onedspec/irsiids/doc/subsets.hlp
@@ -0,0 +1,49 @@
+.help subsets May85 noao.imred.iids/noao.imred.irs
+.ih
+NAME
+subsets - Subtract pairs of spectra in a string
+.ih
+USAGE
+subsets input records
+.ih
+PARAMETERS
+.ls input
+The root file name for the input spectra in the string.
+.le
+.ls records
+The range of spectra indicating the elements of the string.
+The names of the spectra will be formed by appending the range
+elements to the input root name.
+.le
+.ls output
+This is the root file name for the names of the spectra which will
+be created by the subtraction operation.
+.le
+.ls start_rec
+The starting record number to be appended to the root name of the
+created spectra.
+.le
+.ih
+DESCRIPTION
+Pairs of spectra are formed from the input string in the order that
+the record numbers would suggest. 
+The first spectrum in the pair is assumed to be the
+principle spectrum and the second spectrum in the pair is subtracted
+from the first. The result is written out as a new spectrum.
+
+No compensation is made for exposure time during the subtraction.
+The header from the principle spectrum is assigned to the output
+spectrum.
+
+.ih
+EXAMPLES
+The following example forms 50 new spectra from nite1.2001-nite1.2002,
+nite1.2003-nite1.2004, ...
+
+	cl> subsets nite1 2001-2100
+
+The following example creates new spectra from the pairs nite2.2001-nite2.2002,
+nite2.2003-nite2.2004 in spite of the order of the record numbers entered.
+
+	cl> subsets nite2 2001,2003,2002,2004
+.endhelp
diff --git a/noao/onedspec/irsiids/doc/sums.hlp b/noao/onedspec/irsiids/doc/sums.hlp
new file mode 100644
index 00000000..0d8b27e9
--- /dev/null
+++ b/noao/onedspec/irsiids/doc/sums.hlp
@@ -0,0 +1,44 @@
+.help sums Jul85 noao.imred.iids/noao.imred.irs
+.ih
+NAME
+sums -- Generate sums of the sky and object spectra for each aperture
+.ih
+USAGE
+sums input records
+.ih
+PARAMETERS
+.ls input
+The root file name for the input spectra in the string.
+.le
+.ls records
+The range of spectra indicating the elements of the string.
+The names of the spectra will be formed by appending the range
+elements to the input root name.
+.le
+.ls output
+This is the root file name for the names of the spectra which will
+be created by the summation operation.
+.le
+.ls start_rec
+The starting record number to be appended to the root name of the
+created spectra.
+.le
+.ih
+DESCRIPTION
+All the object spectra for each aperture are summed, and the
+sky spectra are also summed to produce two new spectra for
+each observing aperture. Exposure times are accumulated.
+No tests are made to check whether the object is consistent
+among the specified spectra. This could be accomplished by
+checking the titles or telescope positions, but it isn't.
+
+The header parameters OFLAG and BEAM-NUM must be properly
+set in the headers.
+.ih
+EXAMPLES
+The following example forms 4 new spectra from nite1.2001-nite1.2002,
+nite1.2003-nite1.2004, ... assuming this string is derived from
+IIDS spectra.
+
+	cl> sums nite1 2001-2100
+.endhelp
diff --git a/noao/onedspec/irsiids/doc/widstape.hlp b/noao/onedspec/irsiids/doc/widstape.hlp
new file mode 100644
index 00000000..855f223d
--- /dev/null
+++ b/noao/onedspec/irsiids/doc/widstape.hlp
@@ -0,0 +1,90 @@
+.help widstape Mar85 noao.imred.iids/noao.imred.irs
+.ih
+NAME
+widstape -- Write a Cyber style IDSOUT tape
+.ih
+USAGE
+widstape idsout input records
+.ih
+PARAMETERS
+.ls idsout
+The output file name to receive the card-image data. This may be a
+magtape specification (e.g. mta, mtb) or disk file name.
+.le
+.ls input
+The input root file name for the spectra to be written
+.le
+.ls records
+The record string to be appended to the root name to create the image
+names of the spectra to be written.
+.le
+.ls new_tape = no
+If set to yes, the tape is rewound and output begins at BOT. If no,
+output begins at EOT unless an explicit file specification is given
+as part of the magtape file name for parameter "idsout" (e.g. mta[2]).
+If idsout contains a file specification of [1], then writing begins
+at BOT regardless of the value for new_tape.
+.le
+.ls block_size = 3200
+The tape block size in bytes. This must be an integral factor of 80.
+.le
+.ls ebcdic = no
+The default character code is ASCII, but if this parameter is set to yes,
+the output character will be in EBCDIC.
+.le
+.ih
+DESCRIPTION
+The specified spectra are copied to the output file in a card-image format
+defined in the IPPS-IIDS/IRS Reduction Manual. Values from the extended
+image header are used to fill in the observational parameters.
+
+The basic format consists of 4 - 80 byte header cards, 128 data cards
+having 8 data elements per card in 1PE10.3 FORTRAN equivalent format,
+and a trailing blank card for a total of 133 cards. 
+Thus spectra up to 1024 points may be contained in the IDSOUT format. 
+The format is outlined below:
+
+.nf
+ Line	Column	Type
+    1	   1-5	Integer	  Record number within IDSOUT text file
+	  6-10	Integer	  Integration time
+	 11-25	Real	  Wavelength of first bin
+	 26-40	Real	  Dispersion
+	 41-45	Integer	  0 (Index of first pixel)
+	 46-50  Integer	  Line length - 1 (Index of last pixel)
+	 71-80	Integer	  UT time
+    2	  1-10	Real	  Siderial time
+	 11-25	Real	  Right Ascension
+	 26-40	Real	  Declination
+    3	 21-35	Real	  Hour Angle
+	 36-50	Real	  Air mass
+	 51-58	Integer	  UT date
+	 60-76	String	  Image title
+	 78-80	String	  END
+    4	  1-64	String	  Record label
+	 78-80	String	  END
+5-132		Real	  1024 pixel values, 8 per line
+  133			  Blank line
+.fi
+
+The data of type real are in exponent format; i.e FORTRAN 'E' format (1.234e3).
+
+There are no special marks between spectral images, 
+and when multiple spectra are written with a single command, the first card
+of a subsequent spectrum may be within the same physical tape block
+as the last card of the previous spectrum. This assures that all tape
+blocks (except the very last one in the tape file) are all the same
+length.  A double end-of-mark is written after the last spectrum.
+.ih
+EXAMPLES
+The following example writes an IDSOUT format tape starting at the
+beginning of the tape.
+
+	cl> widstape mta nite1 1001-1200 new_tape+
+.ih
+TIME REQUIREMENTS: UNIX/VAX 11/750
+Each spectrum of 1024 points requires about 2 second.
+.ih
+SEE ALSO
+rcardimage, ridsout
+.endhelp
diff --git a/noao/onedspec/irsiids/extinct.cl b/noao/onedspec/irsiids/extinct.cl
new file mode 100644
index 00000000..68c5a2de
--- /dev/null
+++ b/noao/onedspec/irsiids/extinct.cl
@@ -0,0 +1,22 @@
+#{ EXTINCT -- Use the BSWITCH task to perform the correction for
+#            atmospheric extinction.
+
+{
+# Root name
+rt = root
+
+# Records
+rec = records
+
+# Output root
+out = output
+
+# Output starting record
+strt = start_rec
+
+# Do operation
+# Inhibit weighting and statisitic file generation
+
+bswitch  (input=rt, records=rec, output=out, start_rec=strt, subset=nr_aps,
+         weighting=no, ids_mode=no, stats="")
+}
diff --git a/noao/onedspec/irsiids/extinct.par b/noao/onedspec/irsiids/extinct.par
new file mode 100644
index 00000000..a605dae0
--- /dev/null
+++ b/noao/onedspec/irsiids/extinct.par
@@ -0,0 +1,11 @@
+# EXINCT
+
+root,s,a,,,,Root name for spectra file names
+records,s,a,,,,Record string to process
+output,s,a,,,,Root name for spectra to be created
+start_rec,i,a,1,0,9999,Next starting spectral record
+nr_aps,i,h,2
+strt,i,h
+rt,s,h
+out,s,h
+rec,s,h
diff --git a/noao/onedspec/irsiids/flatdiv.par b/noao/onedspec/irsiids/flatdiv.par
new file mode 100644
index 00000000..84de42d4
--- /dev/null
+++ b/noao/onedspec/irsiids/flatdiv.par
@@ -0,0 +1,12 @@
+
+# FLATDIV parameter file
+
+input,s,a,,,,Input image root file name
+records,s,a,,,,Range of spectral records
+output,s,a,,,,Output file root name for new spectra
+start_rec,i,a,1,0,9999,Next starting spectral record
+flat_file,s,a,,,,Image root name for output flat field spectra
+coincor,b,h,)_.coincor,,,Apply coincidence correction to spectra
+ccmode,s,h,)_.ccmode,,,Correction mode (photo|iids)
+deadtime,r,h,)_.deadtime,,,Deadtime in seconds
+power,r,h,)_.power,,,IIDS power law coefficient
diff --git a/noao/onedspec/irsiids/flatfit.par b/noao/onedspec/irsiids/flatfit.par
new file mode 100644
index 00000000..3167697d
--- /dev/null
+++ b/noao/onedspec/irsiids/flatfit.par
@@ -0,0 +1,24 @@
+# FLATFIT parameter file
+
+input,s,a,,,,Input image root file name
+records,s,a,,,,Range of spectral records
+output,s,a,,,,Output file root name for new spectra
+function,s,h,"chebyshev",,,Function to fit (chebyshev|legendre|spline3|spline1)
+order,i,h,6,1,,Fitting order (number of terms)
+niter,i,h,1,1,,Number of rejection iterations
+lower,r,h,2,0,,Lower rejection criterion in sigmas
+upper,r,h,2,0,,Upper rejection criterion in sigmas
+ngrow,i,h,0,,,Growing region
+div_min,r,h,1.0,,,Value to use if division by zero occurs
+interact,b,h,yes,,,Interact with the first accumulation?
+all_interact,b,h,no,,,Interact with all accumulations?
+coincor,b,h,)_.coincor,,,Apply coincidence correction to flats
+ccmode,s,h,)_.ccmode,,,Correction mode (photo|iids)
+deadtime,r,h,)_.deadtime,,,Deadtime in seconds
+power,r,h,)_.power,,,IIDS power law coefficient
+new_order,i,a,4,1,,enter order
+new_lower,r,a,,,,enter nr sigma
+new_upper,r,a,,,,enter nr sigma
+new_niter,i,a,,,,enter nr of iterations
+confirm,b,a,,,,Exit and save solution?
+cursor,*gcur,h,"",,,Graphics cursor input
diff --git a/noao/onedspec/irsiids/getnimage.x b/noao/onedspec/irsiids/getnimage.x
new file mode 100644
index 00000000..c232a85c
--- /dev/null
+++ b/noao/onedspec/irsiids/getnimage.x
@@ -0,0 +1,133 @@
+include	<mach.h>
+
+
+# GET_NEXT_IMAGE -- Use root filename and ranges string (if any) to
+#                  generate the next image filename. Return EOF
+#                  when image list is exhausted.
+
+int procedure get_next_image (infile, records, nrecs, image, sz_name)
+
+int	infile, records[ARB], nrecs, sz_name
+char	image[sz_name]
+
+int	next_num, stat
+int	flag1, flag2, flag3
+char	image_0[SZ_FNAME]
+
+int	clgfil(), get_next_entry(), strlen()
+
+common	/gnicom/ flag1, flag2
+
+data	flag3/YES/
+
+begin
+	# Reset initializer, record counter, and get root name
+	if ((flag1 == YES) || (flag3 == YES)) {
+	    next_num  = -1
+	    call rst_get_entry ()
+	}
+
+	# If no ranges specified, act like template expander
+	if (nrecs == MAX_INT) {
+	    stat = clgfil (infile, image, sz_name)
+
+	# Otherwise append record numbers to first template expansion
+	} else {
+	    if (flag1 == YES) {
+		stat = clgfil (infile, image_0, sz_name)
+		if (stat == EOF)
+		    return (stat)
+	    }
+
+	    stat = get_next_entry (records, next_num)
+	    if (stat != EOF) {
+	 	call strcpy (image_0, image, sz_name)
+		call sprintf (image[strlen(image)+1], sz_name, ".%04d")
+		call pargi (next_num)
+	    }
+	}
+
+	flag1 = NO
+	flag3 = NO
+	return (stat)
+end
+
+
+# Reset the initialization parameter to TRUE
+
+procedure reset_next_image ()
+
+int	flag1, flag2
+common	/gnicom/ flag1, flag2
+
+begin
+	flag1 = YES
+end
+
+
+# GET_NEXT_ENTRY -- Given a list of ranges and the current file number,
+# find and return the next file number in order of entry.  
+# EOF is returned at the end of the list.
+
+int procedure get_next_entry (ranges, number)
+
+int	ranges[ARB]		# Range array
+int	number			# Both input and output parameter
+
+int	ip, first, last, step, next_number, remainder
+int	flag1, flag2, flag3
+
+common	/gnicom/ flag1, flag2
+
+data	flag3/YES/
+
+begin
+	number = number + 1
+	next_number = MAX_INT
+	if ((flag2 == YES) || (flag3 == YES)) {
+	    ip = 1
+	    flag2 = NO
+	    flag3 = NO
+	}
+
+	first = min (ranges[ip], ranges[ip+1])
+	last = max (ranges[ip], ranges[ip+1])
+	step = ranges[ip+2]
+
+	if (number >= first && number <= last) {
+	    remainder = mod (number - first, step)
+	    if (remainder == 0)
+		return (number)
+	    if (number - remainder + step <= last)
+		next_number = number - remainder + step
+	    else
+		go to 10
+
+	} else if (first > number)
+	    next_number = min (next_number, first)
+
+	else {
+10	    ip = ip + 3
+	    if (ranges[ip] != -1 && ranges[ip+1] !=0 && ranges[ip+2] !=0)
+		next_number = min (ranges[ip], ranges[ip+1])
+	}
+
+	if (next_number == MAX_INT) {
+	    ip = 1
+	    flag2 = YES
+	    return (EOF)
+
+	} else {
+	    number = next_number
+	    return (number)
+	}
+end
+
+procedure rst_get_entry ()
+
+int	first, flag2
+common	/gnicom/ first, flag2
+
+begin
+	flag2 = YES
+end
diff --git a/noao/onedspec/irsiids/idsmtn.h b/noao/onedspec/irsiids/idsmtn.h
new file mode 100644
index 00000000..5fa3c73e
--- /dev/null
+++ b/noao/onedspec/irsiids/idsmtn.h
@@ -0,0 +1,101 @@
+# Definitions for the Mountain format IDS tape reader:
+
+define	MAX_RANGES	100
+define	DUMMY		3	# Value returned if DUMMY IDS record is read
+
+define	NBITS_CHAR		(SZB_CHAR * NBITS_BYTE)
+define	SZ_IDS_RECORD		(2108 * 16 / NBITS_CHAR)
+define	NPIX_IDS_REC		1024
+define	LEN_USER_AREA		2880
+define	SZ_IDS_ID		64
+define	NBITS_VN_3WRD_FP	48
+define	NBITS_VN_2WRD_FP	32
+define	NBITS_VN_LONG_INT	32
+define	NBITS_VN_INT		16
+define	NBITS_FORTH_CHAR	8
+define	DATA_BYTE		9			# First byte of data
+define	NBYTES_DATA	(1024 * 32 / NBITS_BYTE)	# Number of data bytes
+define	NBYTES_INT	(NBITS_INT / NBITS_BYTE)
+define	NBYTES_VN_3WRD_FP	6
+define	NBYTES_VN_2WRD_FP	4
+define	NBITS_2WRD_HIGH		8
+define	WRD2_EXP_OFFSET		25
+define	NBITS_2WRD_EXP		6
+define	WRD2_MANT_SIGN		24
+define	WRD2_EXP_SIGN		31
+define	NSIG_VN_BITS		15
+define	VN_LONG_SIGN		31
+define	WRD3_MANT_SIGN		31
+define	WRD3_EXP_SIGN		47
+define	MAX_NCOEFF		25
+
+
+# The control parameter structure is defined below:
+
+define	LEN_CP		10 + SZ_FNAME + 1
+
+define	IS_REDUCED	Memi[$1]
+define	LONG_HEADER	Memi[$1+1]
+define	PRINT_PIXELS	Memi[$1+2]
+define	MAKE_IMAGE	Memi[$1+3]
+define	OFFSET		Memi[$1+4]
+define	DATA_TYPE	Memi[$1+5]
+define	IRAF_FILE	Memc[P2C($1+10)]
+
+
+# The header structure is defined below:
+
+define	LEN_IDS			40 + SZ_IDS_ID + 1
+
+define	HA			Memr[P2R($1)]
+define	AIRMASS			Memr[P2R($1+1)]
+define	RA			Memr[P2R($1+2)]
+define	DEC			Memr[P2R($1+3)]
+define	W0			Memr[P2R($1+4)]
+define	WPC			Memr[P2R($1+5)]
+define	LINE			Memi[$1+6]
+define	NP1			Memi[$1+7]
+define	NP2			Memi[$1+8]
+define	ITM			Memr[P2R($1+9)]
+define	BEAM			Memi[$1+10]
+define	W			Memi[$1+11]
+define	UT			Memr[P2R($1+13)]
+define	ST			Memr[P2R($1+14)]
+define	DF_FLAG			Memi[$1+15]
+define	SM_FLAG			Memi[$1+16]
+define	QF_FLAG			Memi[$1+17]
+define	DC_FLAG			Memi[$1+18]
+define	QD_FLAG			Memi[$1+19]
+define	EX_FLAG			Memi[$1+20]
+define	BS_FLAG			Memi[$1+21]
+define	CA_FLAG			Memi[$1+22]
+define	CO_FLAG			Memi[$1+23]
+define	OFLAG			Memi[$1+24]
+define	POINT			Memi[$1+25]	
+define	DRA			Memi[$1+26]
+define	DDEC			Memi[$1+27]
+define	ALPHA_ID		Memc[P2C($1+35)]
+define	LABEL			Memc[P2C($1+40)]
+
+
+# Bit offsets to various IDS header words are defined below:
+
+define	NREC_OFFSET	((0 * 16) + 1)
+define	RFLAGS_OFFSET	((1 * 16) + 1)
+define	ITM_OFFSET	((2 * 16) + 1)
+define	DATA_OFFSET	((4 * 16) + 1)
+define	W0_OFFSET	((2052 * 16) + 1)
+define	WPC_OFFSET	((2055 * 16) + 1)
+define	NP1_OFFSET	((2058 * 16) + 1)
+define	NP2_OFFSET	((2059 * 16) + 1)
+define	OFLAG_OFFSET	((2060 * 16) + 1)
+define	SMODE_OFFSET	((2061 * 16) + 1)
+define	UT_OFFSET	((2062 * 16) + 1)
+define	ST_OFFSET	((2064 * 16) + 1)
+define	BEAM_OFFSET	((2066 * 16) + 1)
+define	HA_OFFSET	((2067 * 16) + 1)
+define	RA_OFFSET	((2070 * 16) + 1)
+define	DEC_OFFSET	((2073 * 16) + 1)
+define	DRA_OFFSET	((2076 * 16) + 1)
+define	DDEC_OFFSET	((2077 * 16) + 1)
+define	LABEL_OFFSET	((2078 * 16) + 1)
diff --git a/noao/onedspec/irsiids/irsiids.hd b/noao/onedspec/irsiids/irsiids.hd
new file mode 100644
index 00000000..d0a20d98
--- /dev/null
+++ b/noao/onedspec/irsiids/irsiids.hd
@@ -0,0 +1,18 @@
+# Help directory for the IRS/IIDS tasks.
+
+$doc		= "./doc/"
+
+addsets		hlp=doc$addsets.hlp,		src=t_addsets.x
+batchred	hlp=doc$batchred.hlp,		src=batchred.cl
+bswitch		hlp=doc$bswitch.hlp,		src=t_bswitch.x
+coefs		hlp=doc$coefs.hlp,		src=t_coefs.x
+coincor		hlp=doc$coincor.hlp,		src=t_coincor.x
+extinct		hlp=doc$extinct.hlp,		src=extinct.cl
+flatdiv		hlp=doc$flatdiv.hlp,		src=t_flatdiv.x
+flatfit		hlp=doc$flatfit.hlp,		src=t_flatfit.x
+powercor	hlp=doc$powercor.hlp,		src=powercor.cl
+process		hlp=doc$process.hlp,		src=process.cl
+slist1d		hlp=doc$slist1d.hlp,		src=t_slist1d.x
+subsets		hlp=doc$subsets.hlp,		src=t_subsets.x
+sums		hlp=doc$sums.hlp,		src=t_sums.x
+widstape	hlp=doc$widstape.hlp,		src=x_widstape.x
diff --git a/noao/onedspec/irsiids/mkpkg b/noao/onedspec/irsiids/mkpkg
new file mode 100644
index 00000000..01ee3403
--- /dev/null
+++ b/noao/onedspec/irsiids/mkpkg
@@ -0,0 +1,22 @@
+# IRS/IIDS Tasks
+
+$checkout libpkg.a ..
+$update   libpkg.a
+$checkin  libpkg.a ..
+$exit
+
+libpkg.a:
+	coincor.x	
+	conversion.x	
+	getnimage.x	<mach.h>
+	t_addsets.x	<error.h> <imhdr.h>
+	t_bswitch.x	<smw.h> <error.h> <imhdr.h> <mach.h> <time.h>
+	t_coefs.x	<error.h>
+	t_coincor.x	<error.h> <imhdr.h>
+	t_flatdiv.x	<error.h> <imhdr.h>
+	t_flatfit.x	<gset.h> <imhdr.h> <math/curfit.h>
+	t_slist1d.x	<error.h> <fset.h> <imhdr.h> <smw.h>
+	t_subsets.x	<error.h> <imhdr.h>
+	t_sums.x	<error.h> <imhdr.h>
+	t_widstape.x	<error.h> <imhdr.h> <mach.h> <smw.h>
+	;
diff --git a/noao/onedspec/irsiids/powercor.cl b/noao/onedspec/irsiids/powercor.cl
new file mode 100644
index 00000000..a89e2478
--- /dev/null
+++ b/noao/onedspec/irsiids/powercor.cl
@@ -0,0 +1,4 @@
+#{ Apply nonlinear IIDS correction
+
+coincor (input, records, output, start_rec=start_rec, ccmode="iids",
+	power=power)
diff --git a/noao/onedspec/irsiids/powercor.par b/noao/onedspec/irsiids/powercor.par
new file mode 100644
index 00000000..e4e8bb90
--- /dev/null
+++ b/noao/onedspec/irsiids/powercor.par
@@ -0,0 +1,7 @@
+# POWERCOR parameter file
+
+input,s,a,,,,Input image root file name
+records,s,a,,,,Range of spectral records
+output,s,a,,,,Output file root name for new spectra
+start_rec,i,a,1,0,9999,Next starting spectral record
+power,r,h,)_.power,,,Power law coefficient
diff --git a/noao/onedspec/irsiids/slist1d.par b/noao/onedspec/irsiids/slist1d.par
new file mode 100644
index 00000000..ae091a20
--- /dev/null
+++ b/noao/onedspec/irsiids/slist1d.par
@@ -0,0 +1,3 @@
+input,s,a,,,,Input image root file name
+records,s,a,,,,Range of spectral records
+long_header,b,h,no,,,List header in long format
diff --git a/noao/onedspec/irsiids/subsets.par b/noao/onedspec/irsiids/subsets.par
new file mode 100644
index 00000000..5368ca06
--- /dev/null
+++ b/noao/onedspec/irsiids/subsets.par
@@ -0,0 +1,6 @@
+# SUBSETS parameter file
+
+input,s,a,,,,Input image root file name
+records,s,a,,,,Range of spectral records
+output,s,a,,,,Output file root name for new spectra
+start_rec,i,a,1,0,9999,Next starting spectral record
diff --git a/noao/onedspec/irsiids/sums.par b/noao/onedspec/irsiids/sums.par
new file mode 100644
index 00000000..bbeab466
--- /dev/null
+++ b/noao/onedspec/irsiids/sums.par
@@ -0,0 +1,8 @@
+
+# SUMS parameter file
+
+input,s,a,,,,Input image root file name
+records,s,a,,,,Range of spectral records
+output,s,a,,,,Output file root name for new spectra
+start_rec,i,a,1,0,9999,Next starting spectral record
+newoutput,s,q,,,,New output file root name
diff --git a/noao/onedspec/irsiids/t_addsets.x b/noao/onedspec/irsiids/t_addsets.x
new file mode 100644
index 00000000..cd145b4a
--- /dev/null
+++ b/noao/onedspec/irsiids/t_addsets.x
@@ -0,0 +1,195 @@
+include	<error.h>
+include	<imhdr.h>
+
+
+# T_ADDSETS -- Add a series of spectra by subsets. A single spectrum
+#  is produced for every "subset" number of input spectra. The input
+#  list is accumulated until "subset" number of spectra have been
+#  encountered. The result is then written out.
+#
+#  If the input data are calibrated (CA_FLAG = 0) then the result
+#  is an average over the subset size, but the header exposure
+#  time is updated.
+#
+#  If the data is uncalibrated then the resulting spectrum is a sum
+#  of the total counts.
+
+procedure t_addsets ()
+
+pointer	image
+pointer	recstr, ofile
+int	root, start_rec, subset
+int	nrecs
+int	nrem, ifile, ca_flag
+real	itm, expo, wt, wtsum
+bool	weight
+pointer	sp, recs, im, cur_pix, sp_sum
+
+real	imgetr()
+int	clpopni(), clgeti(), imgeti()
+int	get_next_image(), decode_ranges()
+bool	clgetb()
+pointer	immap(), imgl1r()
+
+begin
+	call smark (sp)
+	call salloc (image, SZ_FNAME, TY_CHAR)
+	call salloc (ofile, SZ_FNAME, TY_CHAR)
+	call salloc (recstr, SZ_LINE, TY_CHAR)
+	call salloc (recs, 300, TY_INT)
+
+	# Open input file name template
+	root   = clpopni ("input")
+
+	# Get range specification if any
+	call clgstr ("records", Memc[recstr], SZ_LINE)
+	if (decode_ranges (Memc[recstr], Memi[recs], 100, nrecs) == ERR)
+	    call error (0, "Bad range specification")
+
+	# Get rootname for output files and starting record
+	call clgstr ("output", Memc[ofile], SZ_FNAME)
+	start_rec = clgeti ("start_rec")
+
+	# Get subset size
+	subset = clgeti ("subset")
+
+	# Apply integration time weighting?
+	weight = clgetb ("weighting")
+
+	# Initialize range decoder
+	call reset_next_image ()
+
+	#Initialize file counter
+	ifile = 0
+	wtsum = 0.0
+
+	# Loop over all input images by subsets
+	while (get_next_image (root, Memi[recs], nrecs, Memc[image], 
+	    SZ_FNAME) != EOF) {
+
+	    # Open image
+	    iferr (im = immap (Memc[image], READ_ONLY, 0)) {
+		call erract (EA_WARN)
+		next
+	    }
+
+	    # Allocate space for current subset
+	    if (mod (ifile, subset) == 0) {
+		call calloc (sp_sum, IM_LEN (im,1), TY_REAL)
+
+		# Zero exposure counter
+		expo = 0.0
+	    }
+
+	    # Add in current spectrum
+	    iferr (itm = imgetr (im, "EXPOSURE"))
+		iferr (itm = imgetr (im, "ITIME"))
+		    iferr (itm = imgetr (im, "EXPTIME"))
+			itm = 1
+	    iferr (ca_flag = imgeti (im, "CA-FLAG"))
+		ca_flag = -1
+	    cur_pix = imgl1r (im)
+
+	    # Apply integration time weighting
+	    if (weight)
+		wt = itm
+	    else
+		wt = 1.0
+
+	    if (ca_flag != 0)
+		wt = 1.0
+
+	    wtsum = wtsum + wt
+	    call amulkr (Memr[cur_pix], wt, Memr[cur_pix], IM_LEN(im,1))
+	    call aaddr (Memr[cur_pix], Memr[sp_sum], Memr[sp_sum], IM_LEN(im,1))
+	    expo = expo + itm
+
+	    # Issue status report
+	    call printf ("[%s] added\n")
+		call pargstr (Memc[image])
+
+	    ifile = ifile + 1
+	    if (mod (ifile, subset) == 0) {
+		call wrt_set (Memr[sp_sum], subset, im, Memc[ofile], start_rec,
+		    expo, wtsum, ca_flag)
+		wtsum = 0.0
+		call mfree (sp_sum, TY_REAL)
+	    } else
+		call imunmap (im)
+
+	}
+	# Check that there are no remaining spectra in an unfulfilled subset
+	nrem = mod (ifile, subset)
+	if (nrem != 0) {
+	    call wrt_set (Memr[sp_sum], nrem, im, Memc[ofile], start_rec, 
+		expo, wtsum, ca_flag)
+	    wtsum = 0.0
+	    call mfree (sp_sum, TY_REAL)
+
+	    call eprintf ("Unfulfilled subset accumulation written - ")
+	    call eprintf ("missing %d spectra\n")
+		call pargi (subset - nrem)
+	}
+
+	# Update record number
+	call clputi ("next_rec", start_rec)
+
+	# Free space
+	call sfree (sp)
+	call clpcls (root)
+end
+
+# WRT_SET -- Write spectra ccumulated from the set
+
+procedure wrt_set (sp_sum, subset, im, ofile, start_rec, expo, wtsum, ca_flag)
+
+real	sp_sum[ARB]
+int	subset, start_rec, ca_flag
+pointer	im
+char	ofile[SZ_FNAME]
+real	expo, wtsum
+
+char	newfile[SZ_FNAME]
+pointer imnew, newpix
+
+pointer	impl1r(), immap()
+int	strlen()
+
+begin
+	# Create new spectrum - first make up a name
+	call strcpy (ofile, newfile, SZ_FNAME)
+	call sprintf (newfile[strlen (newfile) + 1], SZ_FNAME, ".%04d")
+	    call pargi (start_rec)
+
+	imnew = immap (newfile, NEW_COPY, im)
+
+	IM_NDIM (imnew) = 1
+	IM_LEN (imnew,1) = IM_LEN (im,1)
+	IM_PIXTYPE (imnew) = TY_REAL
+	call strcpy (IM_TITLE(im), IM_TITLE(imnew), SZ_LINE)
+
+	call imunmap (im)
+
+	newpix = impl1r (imnew)
+
+	# If this spectrum is calibrated, perform an average
+	# weighted by integration time and copy new pixels into image
+	if (ca_flag == 0)
+	    call adivkr (sp_sum, real (wtsum), Memr[newpix], IM_LEN(imnew,1))
+	else
+	    call amovr (sp_sum, Memr[newpix], IM_LEN(imnew,1))
+
+	# Update keyword
+	call imaddr (imnew, "EXPOSURE", expo)
+
+	# Send user report
+	call printf ("writing [%s]: %s\n")
+	    call pargstr (newfile)
+	    call pargstr (IM_TITLE(imnew))
+	call flush (STDOUT)
+
+	call imunmap (imnew)
+
+	# Update record counter
+	start_rec = start_rec + 1
+end
diff --git a/noao/onedspec/irsiids/t_bswitch.x b/noao/onedspec/irsiids/t_bswitch.x
new file mode 100644
index 00000000..0c7b71a5
--- /dev/null
+++ b/noao/onedspec/irsiids/t_bswitch.x
@@ -0,0 +1,924 @@
+include	<error.h>
+include	<imhdr.h>
+include	<mach.h>
+include	<time.h>
+include	<smw.h>
+
+define	MAX_NR_BEAMS	100	# Max number of instrument apertures
+define	MIN_RANGES	100	# Minimum spectra per beam if not given
+
+# T_BSWITCH -- Beam switch a series of spectra to produce a single
+#	     sky subtracted spectrum.
+#
+# The spectra may be extinction corrected if not already done.
+#
+# The summation may include an optional statistical weighting
+# based on the total countrate summed over a user definable
+# piece of the spectrum. If the countrate is <= 0, the
+# spectrum is given zero weight.
+#
+# The data may be organized as data from the IIDS/IRS are usually
+# obtained - where the telescope is beam-switched so that the
+# object is first in one aperture while sky is observed in the other,
+# and then the process is reversed.
+#
+# If the instrument offers many apertures, "nebular" mode can be used
+# to obtain the same effect. Here all apertures observe the object(s)
+# at one time; then the telescope is moved so all apertures are observing
+# sky.
+#
+# Both these methods are considered "idsmode". But if there are a different
+# number of sky observations than object, an imbalance exists. 
+# To account for this possibility, all summations are performed by computing
+# an average countrate over all observations. Sky countrates can then be
+# subtracted from the object. Later the differential countrate is returned
+# to an "equivalent" count by multiplying by the exposure time.
+#
+# Spectra must be dispersion corrected to employ either
+# weighting or extinction correction.
+#
+# The series of spectra may be accumulated in subsets rather than
+# over the entire series by specifying a subset rate. (E.g. for
+# IIDS data a subset rate of 4 would produce a summed pair for
+# every quadruple.)
+
+# Revisions made for WCS support and change from idsmtn.h structure to shdr.h
+# structure.  Because this program is an awful mess the changes were made a
+# small as possible without altering the structure.  (5/1/91, Valdes)
+
+procedure t_bswitch ()
+
+char	image[SZ_FNAME,MAX_NR_BEAMS+1]
+char	rec_numbers[SZ_LINE], title[SZ_LINE,MAX_NR_BEAMS]
+char	ofile[SZ_FNAME], stat_fil[SZ_FNAME]
+int	sfd, nrecsx
+int	i, infile, nrecs, def_beam, start_rec, nimage, sub_rate
+int	records[300], beam_stat[MAX_NR_BEAMS], ncols[MAX_NR_BEAMS]
+bool	idsmode, extinct, stat, weight, eof_test
+pointer	ids[MAX_NR_BEAMS+1] 
+pointer	imnames[MAX_NR_BEAMS]	# Hold pointers to pointers of image names
+pointer	imin, sp, obs
+
+# The following arrays are suffixed by either 'o' for object or 's' for sky
+
+int	ico   [MAX_NR_BEAMS],   ics   [MAX_NR_BEAMS]	# nr obs in beam
+real	expo  [MAX_NR_BEAMS],   exps  [MAX_NR_BEAMS]	# exposure times
+pointer	accumo[MAX_NR_BEAMS+1], accums[MAX_NR_BEAMS+1]	# beam accumulators
+pointer	counto[MAX_NR_BEAMS],   counts[MAX_NR_BEAMS]	# counts in each obs
+
+int	clpopni(), clgeti(), get_next_image(), decode_ranges()
+int	open(), mod()
+pointer	immap()
+bool	clgetb(), streq()
+
+begin
+	call smark (sp)
+	call aclri (ids, MAX_NR_BEAMS+1)
+
+	# Open input filename template
+	infile = clpopni ("input")
+
+	# Get range specification
+	call clgstr ("records", rec_numbers, SZ_LINE)
+	if (decode_ranges (rec_numbers, records, 100, nrecs) == ERR)
+	    call error (0, "Bad range specification")
+
+	# If no ranges is given, filename expansion will occur, so
+	# we must will need some indication of the number of spectra.
+	if (nrecs == MAX_INT)
+	    nrecsx = MIN_RANGES
+	else
+	    nrecsx = nrecs
+
+	# Get root name for new records and starting record number
+	call clgstr ("output", ofile, SZ_FNAME)
+	start_rec = clgeti ("start_rec")
+
+	# Get filename for statistics
+	call clgstr ("stats", stat_fil, SZ_FNAME)
+
+	# Assume spectra are in quadruples?
+	idsmode = clgetb ("ids_mode")
+
+	# Perform de-extinction?
+	extinct = clgetb ("extinct")
+
+	# Use weighting?
+	weight = clgetb ("weighting")
+
+	# Accumulate by subsets? - A very large number implies no subsetting
+	sub_rate = clgeti ("subset")
+
+	# Open statistics file if any
+	if (streq (stat_fil, "")) {
+	    sfd = NULL
+	    stat = false
+	} else {
+	    stat = true
+	    sfd = open (stat_fil, APPEND, TEXT_FILE)
+	}
+
+	# Initialize beam-switch status
+	obs = NULL
+	call init_file (extinct, def_beam, ico, ics, beam_stat)
+
+
+	# Begin cycling through all images - accumulate if possible
+	# by beam number
+
+	# Initialize range decoder
+	call reset_next_image ()
+
+	# Set up for subsets
+	nimage = 0
+	eof_test = false
+
+	repeat {
+
+	while (get_next_image (infile, records, nrecs, image[1,def_beam],
+	    SZ_FNAME) != EOF) {
+
+	    # Attempt to open image with extended header -
+	    iferr (imin = immap (image[1,def_beam], READ_ONLY, 0)) {
+		call eprintf ("[%s]")
+		    call pargstr (image[1,def_beam])
+		call error (0, "Image not found or header info not available")
+	    }
+
+	    # Add in to accumlators
+	    call accum_image (imin, ids, accumo, accums, counto, counts, 
+		ico, ics, expo, exps, image, beam_stat, idsmode, extinct, 
+		weight, nrecsx, ncols, title, imnames, sfd, obs)
+
+	    call printf ("[%s] added\n")
+		call pargstr (image[1,def_beam])
+	    call flush (STDOUT)
+
+	    # Close current image
+	    call imunmap (imin)
+
+	    # Test for subsets
+	    nimage = nimage + 1
+	    if (mod (nimage, sub_rate) == 0)
+		go to 10
+	}
+
+	# Get here by running out of data
+	eof_test = true
+
+	# Must be careful not to write out the last sums if subsets are
+	# in effect because the subset check would already have done so
+	# We can check because "nimage" will not have been bumped
+	# if EOF was encountered.
+
+	if (mod (nimage, sub_rate) != 0) {
+	
+	    # All data has been summed - generate spectra of the accumlations
+10	    call wrt_accum (ids, image, title, accumo, accums, ico, ics, 
+		counto, counts, expo, exps, ncols, beam_stat, idsmode, weight, 
+		extinct, ofile, start_rec, sub_rate)
+
+	    # Generate statistics output for this beam
+	    if (stat)
+		call wrt_stats (sfd, accumo, accums, ico, ics, counto, counts, 
+		    expo, exps, image, beam_stat, title, imnames, weight)
+
+	    # Clear counters and accumulators
+	    call reset_beams (accumo, accums, expo, exps, ico, ics, beam_stat, 
+		ncols)
+	}
+
+	} until (eof_test)
+
+	# Put current record counter back into the parameter file for
+	# subsequent invocations
+	call clputi ("next_rec", start_rec)
+
+	# Close out inputs, outputs, and space
+	do i = 1, MAX_NR_BEAMS+1
+	    call shdr_close (ids[i])
+	if (obs != NULL)
+	    call obsclose (obs)
+	call clpcls (infile)
+	call close (sfd)
+	call sfree (sp)
+end
+
+# ACCUM_IMAGE -- Opens current pixel file, loads header elements,
+#                adds current spectrum to accumulator array(s),
+#                and updates the accumulator status array.
+#		 If not in IDSMODE, then returns both object and
+#		 sky sums for further consideration.
+#		 IDSMODE requires an equal number of each, object and sky,  in
+#		 a sequence of OSSO-OSSO or OSSO-SOOS groups.
+
+procedure accum_image (imin, ids, accumo, accums, counto, counts, ico, ics,
+                expo, exps, image, beam_stat, idsmode, extinct, weight, nrecs,
+		ncols, title, imnames, sfd, obs)
+
+pointer	imin, ids[ARB]
+pointer	imnames[ARB]			# Saved image names for stat printout
+pointer	sfd				# Statistics file
+pointer	obs				# Observatory
+
+pointer	accumo[ARB], accums[ARB]	# Object and sky accumlators
+pointer	counto[ARB], counts[ARB]	#                counting stats
+real	expo  [ARB], exps  [ARB]	#                total exposure times
+int	ico   [ARB], ics   [ARB]	#                number of observations
+
+char	image[SZ_FNAME, MAX_NR_BEAMS+1], title[SZ_LINE,MAX_NR_BEAMS]
+char	observatory[SZ_FNAME]
+int	beam_stat[ARB], ncols[ARB]
+int	dum_beam
+bool	idsmode, extinct, weight, exflag, newobs, obshead
+real	latitude
+
+int	last_len[MAX_NR_BEAMS], name_nr[MAX_NR_BEAMS]
+int	ifile, nr_beams, i, j, def_beam, beam_nr
+int	nwaves, ic, nrecs
+real	airm, wave1, wave2, wt
+pointer	wave_tbl, extn_tbl, ipacc, ipc, mw
+
+real	clgetr(), obsgetr()
+pointer	smw_openim()
+errchk	smw_openim, shdr_open, obsimopen
+
+begin
+	# Bump image file counter
+	ifile = ifile + 1
+
+	# Load header area
+	mw = smw_openim (imin)
+	call shdr_open (imin, mw, 1, 1, INDEFI, SHDATA, ids[def_beam])
+	call smw_close (MW(ids[def_beam]))
+
+	accumo[def_beam] = SY(ids[def_beam])
+
+	# Check for proper flags
+	call flag_chk (ids[def_beam], exflag)
+
+	if (ifile == 1) {
+
+	    # Get region for statistics to operate over -
+	    # Currently only one set of wavelengths is available, but
+	    # at some point, it may be desirable to extend this to
+	    # provide a start and ending wavelength for each aperture
+	    # since an aperture must be considered as an independent
+	    # instrument.
+
+	    # Insert defaults --> entire spectrum
+	    # Now ask user for start and end - if =0.0, use defaults
+	    wave1 = clgetr ("wave1")
+	    wave2 = clgetr ("wave2")
+
+	    if (wave1 == 0.0)
+		wave1 = W0(ids[def_beam])
+	    if (wave2 == 0.0)
+		wave2 = W0(ids[def_beam]) + (IM_LEN(imin,1)-1) * 
+		    WP(ids[def_beam])
+
+	}
+
+	# Determine beam number and add/sub in pixels
+	# Remember that IIDS/IRS "beams" are 0-indexed
+
+	beam_nr = BEAM(ids[def_beam]) + 1
+	if (beam_nr > MAX_NR_BEAMS || beam_nr < 1)
+	    call error (0, "Illegal beam number")
+
+	# Allocate space for this aperture if not already done
+	# Space must be allocated for 2 lines of spectra for
+	# each aperture - Line 1 is used to sum up the most
+	# recent object-sky spectra to maintain the local
+	# statistics. Line 2 is used for the net accumulation
+	# over the entire sequence. The statistics from Line 1
+	# may be used to weigh the observations as they are
+	# added into the Line 2 accumulation.
+	#
+	# For non-IDSMODE the two lines are used for separate
+	# object and sky sums
+
+	if (IS_INDEFI (beam_stat[beam_nr])) {
+	    beam_stat[beam_nr] = 0
+
+	    # Allocate space for the accumulators for this beam nr
+	    call salloc (accumo[beam_nr], IM_LEN(imin,1), TY_REAL)
+	    call salloc (accums[beam_nr], IM_LEN(imin,1), TY_REAL)
+
+	    # Zero object and sky accumulators
+	    call amovkr (0.0, Memr[accumo[beam_nr]], IM_LEN(imin,1))
+	    call amovkr (0.0, Memr[accums[beam_nr]], IM_LEN(imin,1))
+
+
+	    # Allocate space for statistics array - For each beam,
+	    # a series of up to 'nrecs' spectra may be read, and we
+	    # want to keep track of the stats (=countrates) for each
+	    # observation. For non-idsmode, need sky rates too.
+	    call salloc (counto[beam_nr], nrecs, TY_REAL)
+	    if (!idsmode)
+		call salloc (counts[beam_nr], nrecs, TY_REAL)
+
+	    # Allocate space for the image names
+	    call salloc (imnames[beam_nr], nrecs, TY_INT)
+	    name_nr[beam_nr] = 1
+	    do j = 1, nrecs
+		call salloc (Memi[imnames[beam_nr]+j-1], SZ_LINE, TY_CHAR)
+
+	    # Save number of points for checking purposes
+	    last_len[beam_nr] = IM_LEN(imin,1)
+	    ncols[beam_nr] = last_len[beam_nr]
+
+	    # Initialize exposure time
+	    expo[beam_nr] = 0.0
+	    exps[beam_nr] = 0.0
+
+	    nr_beams = nr_beams + 1
+	}
+
+	# If this is an object observation, save the image name
+	if (OFLAG(ids[def_beam]) == 1) {
+	    call strcpy (image[1,def_beam], Memc[Memi[imnames[beam_nr]+
+		name_nr[beam_nr]-1]], SZ_LINE)
+	    name_nr[beam_nr] = name_nr[beam_nr] + 1
+	}
+
+	# If an object observation, save the header elements --
+	# NOTE that if we get >1 objects before getting a sky, only
+	# the last  observation header is saved!
+
+	# The pixel data will be the sum of all objects until the
+	# |object-sky| count = 0 -- Thus, beam switching does not
+	# necessarily accumulate by pairs, but depends on how the
+	# sequence of observations are presented to the program.
+
+	# The following test has been deleted so that headers
+	# will be saved for sky frames as well. This is necessary
+	# if BSWITCH is to perform the function of EXTINCTION
+	# only when sky frames are to be written out as well.
+
+	if (OFLAG(ids[def_beam]) == 1 || !idsmode) {
+	    # Save headers - could probably be done faster by AMOV
+	    call shdr_copy (ids[def_beam], ids[beam_nr], NO)
+
+	    # Fix airmass if necessary
+	    if (extinct && IS_INDEF (AM(ids[beam_nr]))) {
+		call clgstr ("observatory", observatory, SZ_FNAME)
+		call obsimopen (obs, imin, observatory, NO, newobs, obshead)
+		if (newobs) {
+		    call obslog (obs, "BSWITCH", "latitude", STDOUT)
+		    if (sfd != NULL)
+			call obslog (obs, "BSWITCH", "latitude", sfd)
+		}
+		latitude = obsgetr (obs, "latitude")
+		call get_airm (RA(ids[beam_nr]), DEC(ids[beam_nr]),
+		    HA(ids[beam_nr]), ST(ids[beam_nr]), latitude,
+		    AM(ids[beam_nr]))
+	    }
+
+	    call strcpy (image[1,def_beam], image[1,beam_nr], SZ_FNAME)
+
+	    # Save length - Each beam may be independent sizes
+	    ncols[beam_nr] = IM_LEN(imin,1)
+
+	    # Save title, too for same reason
+	    call strcpy (IM_TITLE(imin), title[1,beam_nr], SZ_LINE)
+	}
+
+	# Verify length
+	if (last_len[beam_nr] != ncols[beam_nr]) {
+	    call eprintf ("[%s] -- Length not consistent %d\n")
+		call pargstr (image[1,beam_nr])
+		call pargi   (ncols[beam_nr])
+	    ncols[beam_nr] = min (ncols[beam_nr], last_len[beam_nr])
+	}
+	last_len[beam_nr] = ncols[beam_nr]
+
+
+	# Check to see if a pair is obtained - then perform statistics
+	# and add into global accumulator
+
+	if (idsmode) {
+
+	    # Add spectrum to local accumulation buffer --> Use SKY buffer
+	    # At this point of deriving a sequentially local sum, weighting
+	    # is not used. 
+
+	    call add_spec (Memr[accumo[def_beam]], Memr[accums[beam_nr]],
+		beam_stat[beam_nr], OFLAG(ids[def_beam]), last_len[beam_nr])
+
+	    # IDSMODE requires that every 2N observations produce an
+	    # OBJECT-SKY pair
+	    if (mod (ifile, 2*nr_beams) == 0)
+
+		# Review all beams in use for non-zero pairings
+		do i = 1, MAX_NR_BEAMS
+		    if (!IS_INDEFI (beam_stat[i]) && beam_stat[i] != 0)
+			call error (0, "Spectra are not in quadruples")
+
+
+	    # Object and sky exposure times must be equal.
+	    if (OFLAG(ids[def_beam]) == 1) {
+		expo[beam_nr] = expo[beam_nr] + IT(ids[def_beam])
+
+		# Increment number of object observations for this beam
+		ico[beam_nr] = ico[beam_nr] + 1
+	    }
+
+
+	    if (beam_stat[beam_nr] == 0) {
+		# Add up all counts within a region for statistics of objects
+		# This must be kept separately for each beam number and for
+		# each observation
+
+		# First convert to counts per second (CPS)
+		call adivkr (Memr[accums[beam_nr]], IT(ids[def_beam]),
+		    Memr[accums[beam_nr]], last_len[beam_nr])
+
+		# Sum CPS in statistics region
+		call sum_spec (Memr[accums[beam_nr]], wave1, wave2, 
+		    W0(ids[def_beam]), WP(ids[def_beam]), Memr[counto[beam_nr]+
+		    ico[beam_nr]-1], last_len[beam_nr])
+
+		# De-extinct spectrum
+		if (extinct && !exflag) {
+		    airm = AM(ids[beam_nr])
+		    call de_ext_spec (Memr[accums[beam_nr]], airm, 
+			W0(ids[def_beam]), WP(ids[def_beam]), Memr[wave_tbl], 
+			Memr[extn_tbl], nwaves, last_len[beam_nr])
+		}
+
+		# Add to global accumulator
+		# Use weights which are proportional to countrate, if desired
+		if (weight) {
+		    wt = Memr[counto[beam_nr]+ico[beam_nr]-1]
+		    call amulkr (Memr[accums[beam_nr]], wt,
+			Memr[accums[beam_nr]], last_len[beam_nr])
+		}
+
+		# And add into global sum
+		call aaddr (Memr[accums[beam_nr]], Memr[accumo[beam_nr]],
+		    Memr[accumo[beam_nr]], last_len[beam_nr])
+	    }
+
+	} else {
+	    # Non IDSMODE -accumulate separate object and sky CPS sums
+
+	    # Set pointers and update obj-sky parameters
+	    if (OFLAG(ids[def_beam]) == 1) {
+		beam_stat[beam_nr] = beam_stat[beam_nr] + 1
+		ipacc = accumo[beam_nr]
+		ipc   = counto[beam_nr]
+		ico[beam_nr] = ico[beam_nr] + 1
+		ic = ico[beam_nr]
+		expo[beam_nr] = expo[beam_nr] + IT(ids[def_beam])
+	    } else {
+		beam_stat[beam_nr] = beam_stat[beam_nr] - 1
+		ipacc = accums[beam_nr]
+		ipc   = counts[beam_nr]
+		ics[beam_nr] = ics[beam_nr] + 1
+		ic = ics[beam_nr]
+		exps[beam_nr] = exps[beam_nr] + IT(ids[def_beam])
+	    }
+
+	    # First convert to counts per second (CPS)
+	    call adivkr (Memr[accumo[def_beam]], IT(ids[def_beam]),
+		    Memr[accumo[def_beam]], last_len[beam_nr])
+
+	    # Get counting stats
+	    call sum_spec (Memr[accumo[def_beam]], wave1, wave2, 
+		W0(ids[def_beam]), WP(ids[def_beam]), Memr[ipc+ic-1], 
+		last_len[beam_nr])
+
+	    # De-extinct spectrum
+	    if (extinct && !exflag) {
+		airm = AM(ids[beam_nr])
+		call de_ext_spec (Memr[accumo[def_beam]], airm, 
+		    W0(ids[def_beam]), WP(ids[def_beam]), Memr[wave_tbl], 
+		    Memr[extn_tbl], nwaves, last_len[beam_nr])
+	    }
+
+	    if (weight) {
+		wt = Memr[ipc+ic-1]
+		call amulkr (Memr[accumo[def_beam]], wt, Memr[accumo[def_beam]],
+		    last_len[beam_nr])
+	    }
+
+	    # Add into appropriate accumulator
+	    call aaddr (Memr[accumo[def_beam]], Memr[ipacc], Memr[ipacc],
+		last_len[beam_nr])
+	}
+	
+	return
+
+# INIT_FILE -- Zero the file initializer, the beam counter, beam stats
+#	       and read the extinction data if necessary
+
+entry init_file (extinct, dum_beam, ico, ics, beam_stat)
+
+	ifile = 0
+	nr_beams = 0
+	def_beam = MAX_NR_BEAMS + 1
+	dum_beam = def_beam
+
+	do i = 1, MAX_NR_BEAMS {
+	    beam_stat[i] = INDEFI
+	    ico[i] = 0
+	    ics[i] = 0
+	}
+
+	# If extinction required, read in extinction file, and sensitivity file
+	if (extinct)
+	    call get_extn (wave_tbl, extn_tbl, nwaves)
+
+	return
+
+# INIT_NAME -- Reset name index counter for a beam number
+
+entry init_name (dum_beam)
+
+	name_nr[dum_beam] = 1
+	return
+end
+
+# ACCUM_OUT -- Checks accumulator flags and writes out a new summed
+#              image if the count is zero
+
+procedure accum_out (accum, image, ncols, title, root, rec, beam_nr,
+	bsflag, itm, exflag)
+
+real	accum[ARB], itm
+char	image[SZ_FNAME], title[SZ_LINE], root[SZ_FNAME]
+int	ncols, rec, beam_nr
+int	bsflag, exflag
+
+pointer	imin, imout, spec
+char	bs_image[SZ_FNAME]
+
+pointer	immap(), impl1r()
+
+begin
+	# Create new image with user area
+	# Use ROOT for spectrum name and increment starting record number
+
+	call sprintf (bs_image, SZ_FNAME, "%s.%04d")
+	    call pargstr (root)
+	    call pargi (rec)
+
+	rec = rec + 1
+
+	# Provide user info
+	call printf ("writing: [%s] %s\n")
+	    call pargstr (bs_image)
+	    call pargstr (title)
+	call flush (STDOUT)
+
+	imin = immap (image, READ_ONLY, 0)
+	imout = immap (bs_image, NEW_COPY, imin)
+
+	# Add standard image header
+	IM_NDIM(imout) = 1
+	IM_LEN(imout,1) = ncols
+	IM_PIXTYPE(imout) = TY_REAL
+	call strcpy (title, IM_TITLE(imout), SZ_LINE)
+
+	# Write out pixels
+	spec = impl1r (imout)
+	call amovr (accum, Memr[spec], ncols)
+
+	# Update changed parameters
+	if(bsflag == 1)
+	    call imaddi (imout, "BS-FLAG", bsflag)
+	call imaddr (imout, "EXPTIME", itm)
+	call imaddi (imout, "EX-FLAG", exflag)
+
+	call imunmap (imin)
+	call imunmap (imout)
+
+	# Store new image name back into image
+	call strcpy (bs_image, image, SZ_FNAME)
+end
+
+# ACCUM_NORM - Normalize weighted rate and convert to counts
+
+procedure accum_norm (accum, nr, counts, exp, ncols, weight)
+
+real	accum[ARB], counts[ARB], exp
+int	nr, ncols
+bool	weight
+
+real	sum_wt
+int	i
+
+begin
+	# The accumulation is an array weighted by non-normalized weights
+	# Normalize to total weight to produce a true weighted average
+	# and multiply by the total exposure to produce
+	# an equivalent sum
+
+	# Add up all weighting factors
+	if (weight) {
+	    sum_wt = 0.0
+	    do i = 1, nr
+		sum_wt = sum_wt + counts[i]
+	} else
+	    sum_wt = real (nr)
+
+	if (sum_wt == 0.0)
+	    sum_wt = 1.0
+
+	# Correct for exposure time
+	sum_wt = exp / sum_wt
+
+	call amulkr (accum, sum_wt, accum, ncols)
+end
+
+# WRT_ACCUM -- Write out accumulations as spectra
+
+procedure wrt_accum (ids, image, title, accumo, accums, ico, ics, 
+	counto, counts, expo, exps, ncols, beam_stat, idsmode, weight, 
+	extinct, ofile, start_rec, sub_rate)
+
+pointer	ids[ARB]
+char	image[SZ_FNAME,MAX_NR_BEAMS+1], title[SZ_LINE,MAX_NR_BEAMS]
+
+pointer	accumo[ARB], accums[ARB]
+pointer	counto[ARB], counts[ARB]
+int	ico   [ARB], ics   [ARB]
+real	expo  [ARB], exps  [ARB]
+int	ncols[ARB]
+int	beam_stat[ARB]
+bool	idsmode, weight, extinct
+char	ofile[SZ_FNAME]
+int	start_rec, sub_rate, bsflag
+
+int	i, nr_beams
+real	exp_ratio
+
+begin
+	# First compute number of beams
+	nr_beams = 0
+	do i = 1, MAX_NR_BEAMS
+	    if (!IS_INDEFI (beam_stat[i]) && ((ico[i] > 0) || (ics[i] > 0)))
+		nr_beams = nr_beams + 1
+
+	# For all present apertures, write out a spectrum
+	do i = 1, MAX_NR_BEAMS {
+
+	    if (!IS_INDEFI (beam_stat[i]) && ((ico[i] > 0) || (ics[i] > 0))) {
+		if (beam_stat[i] != 0 && idsmode) {
+		    call eprintf ("Non-equal number of obj-sky observations")
+		    call eprintf (" beam: %d - residual: %d\n")
+			call pargi (i-1)
+		        call pargi (beam_stat[i])
+
+		    # Reset to 0 and force output
+		    beam_stat[i] = 0
+		}
+
+		# The accumulator has a total CPS using non-normalized
+		# weights - apply normalization and exposure time to
+		# generate an equivalent COUNT sum.
+		call accum_norm (Memr[accumo[i]], ico[i], Memr[counto[i]], 
+			expo[i], ncols[i], weight)
+
+		if (!idsmode) {
+		    # Separate object and sky sums require sky info
+		    call accum_norm (Memr[accums[i]], ics[i], Memr[counts[i]],
+			exps[i], ncols[i], weight)
+
+		    # Then normalize sky exposure time to that of object
+		    if (exps[i] != 0.0)
+			exp_ratio = expo[i]/exps[i]
+		    else
+			exp_ratio = 1.0
+
+		    # Check that some object observtion was made
+		    # If not, then we only have sky data so multiply by -1
+		    # so that the subsequent subtraction will produce a
+		    # positive sky
+		    if (expo[i] == 0.0)
+			exp_ratio = -1.0
+
+		    if (exp_ratio != 1.0)
+			call amulkr (Memr[accums[i]], exp_ratio,
+			    Memr[accums[i]], ncols[i])
+
+		    # Finally subtract sky from object equivalent counts
+		    call asubr (Memr[accumo[i]], Memr[accums[i]], 
+			Memr[accumo[i]], ncols[i])
+
+		}
+		# Set header flags
+		# BS flag is not set if the subset rate equals the
+		# number of apertures since each record in is copied out
+		if (sub_rate > nr_beams)
+		    bsflag = 1
+		else
+		    bsflag = -1
+
+		if (OFLAG(ids[i]) == 1)
+		    IT(ids[i]) = expo[i]
+		else
+		    IT(ids[i]) = exps[i]
+
+		if (extinct)
+		    EC(ids[i]) = 0
+
+		# And write out spectrum, at last
+		call accum_out (Memr[accumo[i]], image[1,i], 
+		    ncols[i], title[1,i], ofile, start_rec, i,
+		    bsflag, IT(ids[i]), EC(ids[i]))
+
+		# Reset name entry counter
+		call init_name (i)
+	    }
+
+	}
+end
+
+# RESET_BEAMS -- Zeroes the counters and accumulators for additional
+#                cases
+
+procedure reset_beams (accumo, accums, expo, exps, ico, ics, beam_stat, ncols)
+
+pointer	accumo[ARB], accums[ARB]
+real	expo  [ARB], exps  [ARB]
+int	ico   [ARB], ics   [ARB]
+int	beam_stat[ARB], ncols[ARB]
+
+int	i
+
+begin
+	do i = 1, MAX_NR_BEAMS
+	    if (!IS_INDEFI (beam_stat[i])) {
+
+		expo[i] = 0.0
+		exps[i] = 0.0
+		ico[i]  = 0
+		ics[i]  = 0
+		call amovkr (0.0, Memr[accumo[i]], ncols[i])
+		call amovkr (0.0, Memr[accums[i]], ncols[i])
+	    }
+end
+
+# WRT_STATS -- Write out statistics file
+
+procedure wrt_stats (fd, accumo, accums, ico, ics, counto, counts,
+     expo, exps, image, beam_stat, title, imnames, weight)
+
+int	fd
+pointer	accumo[ARB], accums[ARB], counto[ARB], counts[ARB]
+real	expo[ARB], exps[ARB]
+int	ico[ARB], ics[ARB], beam_stat[ARB]
+char	image[SZ_FNAME,MAX_NR_BEAMS+1]
+char	title[ARB]
+pointer	imnames[ARB]
+bool	weight
+
+int	i, j
+real	cmaxo, cmaxs
+char	ctime[SZ_TIME]
+
+long	clktime()
+
+begin
+	# Issue time stamp
+	call cnvtime (clktime (long(0)), ctime, SZ_TIME)
+	call fprintf (fd, "%s\n\n")
+	    call pargstr (ctime)
+
+	# Issue message if weighted sums in effect
+	if (weight)
+	    call fprintf (fd, "--> Using weighted averages <--\n\n")
+
+	# Cycle over beams
+	do i = 1, MAX_NR_BEAMS {
+	    if (!IS_INDEFI (beam_stat[i])) {
+
+		# Write out Object stats if any
+		if (ico[i] > 0) {
+		    call fprintf (fd, "Object statistics for beam %d -->[%s]\n")
+			call pargi (i-1)
+			call pargstr (image[1,i])
+		    call fprintf (fd, "Title: %s\n")
+			call pargstr (title)
+
+		    # Find maximum count value for this beam
+		    cmaxo = Memr[counto[i]]
+
+		    do j = 1, ico[i]
+		        cmaxo = max (cmaxo, Memr[counto[i]+j-1])
+
+		    call fprintf (fd,
+			"Obs  Relative CPS   Image%12wPeak CPS = %10.3g\n")
+			call pargr (cmaxo)
+
+		    if (cmaxo == 0.0)
+			cmaxo = 1.0
+
+		    do j = 1, ico[i] {
+			call fprintf (fd, "%3d    %8.2f    [%s]\n")
+			    call pargi (j)
+			    call pargr (Memr[counto[i]+j-1] / cmaxo)
+			    call pargstr (Memc[Memi[imnames[i]+j-1]])
+		    }
+		}
+
+		# Write out sky stats if any
+		if (ics[i] > 0) {
+		    call fprintf (fd, "Sky statistics for beam %d\n")
+			call pargi (i-1)
+
+		    cmaxs = Memr[counts[i]]
+
+		    do j = 1, ics[i]
+		        cmaxs = max (cmaxs, Memr[counts[i]+j-1])
+
+		    call fprintf (fd, "Obs  Relative CPS   Peak CPS = %10.3g\n")
+			call pargr (cmaxs)
+
+		    if (cmaxs == 0.0)
+			cmaxs = 1.0
+
+		    do j = 1, ics[i] {
+			call fprintf (fd, "%3d    %8.2f\n")
+			    call pargi (j)
+			    call pargr (Memr[counts[i]+j-1] / cmaxs)
+		    }
+		}
+
+		call fprintf (fd, "\n\n")
+	    }
+	}
+end
+
+
+# ADD_SPEC -- Accumulate spectrum into array - either add or subtract
+#             Returns status = net number of object - sky apectra
+#                            = 0 for equal numbers to indicate further
+#                              processing may take place
+
+procedure add_spec (inspec, accum, stat, flag, len)
+
+real	inspec[ARB], accum[ARB]
+int	stat, flag, len
+
+int	i, add_sub
+
+begin
+	add_sub = 0
+
+	# Is this an Object or Sky?
+	# If flag is neither 0 or 1, spectrum is ignored
+	if (flag == 0)
+	    add_sub = -1
+	if (flag == 1)
+	    add_sub = +1
+
+	if (add_sub == 0) {
+	    stat = INDEFI
+	    return
+	}
+
+	# Is accumulator to be cleared?
+	if (IS_INDEFI (stat) || stat == 0) {
+	    call amulkr (inspec, real (add_sub), accum, len)
+	    stat = add_sub
+
+	} else {
+	    # Add into accumulator
+	    do i = 1, len
+		accum[i] = accum[i] + add_sub * inspec[i]
+
+	    stat = stat + add_sub
+	}
+end
+
+# FLAG_CHK -- Check header flags prior to beam switching
+
+procedure flag_chk (ids, exflag)
+
+pointer	ids
+bool	exflag
+
+int	bsflag, imgeti()
+
+begin
+	# BS requires 
+	# 1. dispersion corrected spectra
+	# 2. non-beam switched
+	# 3. may be either extinction corrected or not
+
+	if (DC(ids) != DCLINEAR)
+	    call error (0, "Spectrum not dispersion corrected")
+
+	iferr (bsflag = imgeti (IM(ids), "BS-FLAG"))
+	    bsflag = -1
+	if (bsflag == 1)
+	    call error (0, "Spectrum already beam-switched")
+
+	if (EC(ids) == ECYES)
+	    exflag = true
+	else
+	    exflag = false
+end
diff --git a/noao/onedspec/irsiids/t_coefs.x b/noao/onedspec/irsiids/t_coefs.x
new file mode 100644
index 00000000..656e777d
--- /dev/null
+++ b/noao/onedspec/irsiids/t_coefs.x
@@ -0,0 +1,88 @@
+include	<error.h>
+
+# COEFS --  Convert IIDS/IRS coeffients to IDENTIFY database entry.
+
+procedure t_coefs ()
+
+int	root			# List of input root names
+pointer	database		# Output database directory
+
+int	i, nrecs, ncoefs
+real	coef
+pointer	sp, image, dtname, recs, im, dt
+
+real	imgetr()
+int	clpopni(), imgeti(), get_next_image(), decode_ranges()
+pointer	immap(), dtmap1()
+errchk	imgetr, dtmap1
+
+begin
+	call smark (sp)
+	call salloc (image, SZ_LINE, TY_CHAR)
+	call salloc (database, SZ_FNAME, TY_CHAR)
+	call salloc (dtname, SZ_FNAME, TY_CHAR)
+	call salloc (recs, 300, TY_INT)
+
+	root   = clpopni ("input")
+	call clgstr ("records", Memc[image], SZ_LINE)
+	call clgstr ("database", Memc[database], SZ_LINE)
+
+	if (decode_ranges (Memc[image], Memi[recs], 100, nrecs) == ERR)
+	    call error (0, "Bad range specification")
+
+	# Loop over all input images - print name on STDOUT
+	call reset_next_image ()
+	while (get_next_image (root, Memi[recs], nrecs, Memc[image],
+	    SZ_LINE) != EOF) {
+	    iferr (im = immap (Memc[image], READ_ONLY, 0)) {
+		call erract (EA_WARN)
+		next
+	    }
+
+	    iferr (ncoefs = imgeti (im, "DF-FLAG"))
+		ncoefs = -1
+	    if (ncoefs > 1) {
+		call strcpy ("id", Memc[dtname], SZ_FNAME)
+		call imgcluster (Memc[image], Memc[dtname+2], SZ_FNAME)
+		dt = dtmap1 (Memc[database], Memc[dtname], APPEND)
+
+		call dtptime (dt)
+		call dtput (dt, "begin\tidentify %s\n")
+		    call pargstr (Memc[image])
+		call dtput (dt, "\tid\t%s\n")
+		    call pargstr (Memc[image])
+		call dtput (dt, "\ttask\tidentify\n")
+		call dtput (dt, "\timage\t%s\n")
+		    call pargstr (Memc[image])
+
+		# Convert coefficients
+		call dtput (dt, "\tcoefficients\t%d\n")
+		    call pargi (ncoefs+4)
+		call dtput (dt, "\t\t2\n")
+		call dtput (dt, "\t\t%1d\n")
+		    call pargi (ncoefs)
+		call dtput (dt, "\t\t1\n")
+		call dtput (dt, "\t\t1024\n")
+
+		do i = 1, ncoefs {
+		    call sprintf (Memc[dtname], SZ_FNAME, "DF%d")
+			call pargi (i)
+		    coef = imgetr (im, Memc[dtname])
+		    call dtput (dt, "\t\t%10.4f\n")
+			call pargr (coef)
+		}
+
+		call dtput (dt, "\n")
+		call dtunmap (dt)
+	    }
+
+	    call printf ("[%s] %d coefficients written\n")
+		call pargstr (Memc[image])
+		call pargi (max (0, ncoefs))
+	    call flush (STDOUT)
+	    call imunmap (im)
+	}
+
+	call clpcls (root)
+	call sfree (sp)
+end
diff --git a/noao/onedspec/irsiids/t_coincor.x b/noao/onedspec/irsiids/t_coincor.x
new file mode 100644
index 00000000..ad2d8bd4
--- /dev/null
+++ b/noao/onedspec/irsiids/t_coincor.x
@@ -0,0 +1,102 @@
+include	<error.h>
+include	<imhdr.h>
+
+
+# T_COINCOR -- Apply coincidence corrections to spectra
+
+procedure t_coincor ()
+
+int	root, start_rec, ccmode, npts, nrecs, coflag
+real	dtime, power, expo
+pointer	sp, image, ofile, str, recs, imin, imout, pixin, pixout
+
+int	clpopni(), clgeti(), clgwrd(), imgeti()
+int	get_next_image(), decode_ranges()
+real	clgetr(), imgetr()
+pointer	immap(), imgl1r(), impl1r()
+errchk	coincor
+
+begin
+	# Allocate memory
+	call smark (sp)
+	call salloc (image, SZ_FNAME, TY_CHAR)
+	call salloc (ofile, SZ_FNAME, TY_CHAR)
+	call salloc (str, SZ_LINE, TY_CHAR)
+	call salloc (recs, 300, TY_INT)
+
+	# Get parameters
+	root   = clpopni ("input")
+	call clgstr ("output", Memc[ofile], SZ_FNAME)
+	if (Memc[ofile] != EOS)
+	    start_rec = clgeti ("start_rec")
+	ccmode = clgwrd ("ccmode", Memc[str], SZ_LINE, ",photo,iids,")
+	dtime = clgetr ("deadtime")
+	power = clgetr ("power")
+	call clgstr ("records", Memc[str], SZ_LINE)
+
+	# Initialize
+	if (decode_ranges (Memc[str], Memi[recs], 100, nrecs) == ERR)
+	    call error (0, "Bad range specification")
+	call reset_next_image ()
+
+	# Loop over all input images by subsets
+	while (get_next_image (root, Memi[recs], nrecs, Memc[image],
+	    SZ_FNAME) != EOF) {
+
+	    # Open input image and check coincidence flag
+	    iferr (imin = immap (Memc[image], READ_WRITE, 0)) {
+		call erract (EA_WARN)
+	        start_rec = start_rec + 1
+		next
+	    }
+	    iferr (coflag = imgeti (imin, "CO-FLAG"))
+		coflag = -1
+	    if (coflag > 0) {
+		call printf ("[%s] already coincidence corrected\n")
+		    call pargstr (IM_HDRFILE(imin))
+		call flush (STDOUT)
+		call imunmap (imin)
+		next
+	    }
+
+	    # Open output image
+	    if (Memc[ofile] != EOS) {
+		call sprintf (Memc[str], SZ_LINE, "%s.%04d")
+		    call pargstr (Memc[ofile])
+		    call pargi (start_rec)
+		start_rec = start_rec + 1
+
+		imout = immap (Memc[str], NEW_COPY, imin)
+		IM_PIXTYPE (imout) = TY_REAL
+	    } else
+		imout = imin
+
+	    # Apply coincidence correction
+	    pixin = imgl1r (imin)
+	    pixout = impl1r (imout)
+	    npts = IM_LEN (imin, 1)
+	    iferr (expo = imgetr (imin, "EXPOSURE"))
+		iferr (expo = imgetr (imin, "ITIME"))
+		    iferr (expo = imgetr (imin, "EXPTIME"))
+			expo = 1
+	    call coincor (Memr[pixin], Memr[pixout], npts, expo, coflag,
+		dtime, power, ccmode)
+
+	    # Update flag and write status
+	    call imaddi (imout, "CO-FLAG", coflag)
+	    call printf ("[%s] --> [%s] %s\n")
+		call pargstr (IM_HDRFILE(imin))
+		call pargstr (IM_HDRFILE(imout))
+		call pargstr (IM_TITLE(imout))
+	    call flush (STDOUT)
+
+	    # Close images
+	    if (imout != imin)
+		call imunmap (imout)
+	    call imunmap (imin)
+	}
+
+	call clputi ("next_rec", start_rec)
+	call clpcls (root)
+	call sfree (sp)
+end
diff --git a/noao/onedspec/irsiids/t_flatdiv.x b/noao/onedspec/irsiids/t_flatdiv.x
new file mode 100644
index 00000000..37186878
--- /dev/null
+++ b/noao/onedspec/irsiids/t_flatdiv.x
@@ -0,0 +1,276 @@
+include	<imhdr.h>
+include	<error.h>
+
+define	MAX_NR_BEAMS	100	# Max number of instrument apertures
+
+# T_FLATDIV -- Divide by a flat field spectrum. This is basically
+#  a simple division of two vectors but with the following
+#  additional functions:
+#
+#	1. Check the processing flag of the input spectra to avoid
+#	   double processing, and set the flag if the processing is
+#	   performed.
+#	2. Trap division by zero errors
+#	3. Optionally apply coincidence corrections
+
+procedure t_flatdiv ()
+
+int	root, start_rec
+int	nrecs
+int	len_flat
+int	ccmode, qd_flag
+real	dtime
+real	power
+bool	coincidence
+pointer	sp, image, str, ofile, flat, recs, bstat, flatsp, im
+
+int	clpopni(), clgeti(), clgwrd(), imgeti()
+int	get_next_image(), decode_ranges()
+real	clgetr()
+bool	clgetb()
+pointer	immap()
+errchk	get_flatsp
+
+begin
+	call smark (sp)
+	call salloc (image, SZ_FNAME, TY_CHAR)
+	call salloc (str, SZ_LINE, TY_CHAR)
+	call salloc (ofile, SZ_FNAME, TY_CHAR)
+	call salloc (flat, SZ_FNAME, TY_CHAR)
+	call salloc (recs, 300, TY_INT)
+	call salloc (bstat, MAX_NR_BEAMS, TY_INT)
+
+	# Open input file name template
+	root   = clpopni ("input")
+
+	# Get range specification if any
+	call clgstr ("records", Memc[str], SZ_LINE)
+	if (decode_ranges (Memc[str], Memi[recs], 100, nrecs) == ERR)
+	    call error (0, "Bad range specification")
+
+	# Get rootname for output files and starting record
+	# Subtract 1 from start_rec because 1 will be added later.
+	call clgstr ("output", Memc[ofile], SZ_FNAME)
+	start_rec = clgeti ("start_rec") - 1
+
+	# Get flat field spectrum root name
+	call clgstr ("flat_file", Memc[flat], SZ_FNAME)
+
+	# Apply coincidence corrections?
+	coincidence = clgetb ("coincor")
+	if (coincidence) {
+	    ccmode = clgwrd ("ccmode", Memc[str], SZ_LINE, ",photo,iids,")
+	    dtime  = clgetr ("deadtime")
+	    power = clgetr ("power")
+	}
+
+	# Initialize beam number status
+	call init_bs (Memi[bstat])
+
+	# Initialize range decoder
+	call reset_next_image ()
+
+	# Loop over all input images - divide and make new image.
+	# The output record number is incremented in all cases.
+	while (get_next_image (root, Memi[recs], nrecs, Memc[image],
+	    SZ_FNAME) != EOF) {
+	    start_rec = start_rec + 1
+
+	    # Open image
+	    iferr (im = immap (Memc[image], READ_ONLY, 0)) {
+		call erract (EA_WARN)
+		next
+	    }
+
+	    # Get header
+	    iferr (qd_flag = imgeti (im, "QD-FLAG"))
+		qd_flag = -1
+	    
+	    # Verify divide flag
+	    if (qd_flag != 0) {
+
+		# Get flat field spectrum if needed
+		call get_flatsp (im, flatsp, Memc[flat], Memi[bstat], len_flat)
+
+		# Calibrate the current spectrum and make a calibrated version
+		call divide (im, flatsp, len_flat, Memc[image], Memc[ofile],
+		    start_rec, coincidence, ccmode, dtime, power)
+
+	    } else {
+		call eprintf ("[%s] already divided - ignored\n")
+		    call pargstr (Memc[image])
+	    }
+	}
+
+	# Update record number
+	call clputi ("next_rec", start_rec)
+
+	# Free space
+	call sfree (sp)
+	call clpcls (root)
+end
+
+# GET_FLATSP -- Load flat field spectrum for the current beam number
+
+procedure get_flatsp (im, sp, flat_file, beam_stat, len_flat)
+
+pointer	im, sp
+char	flat_file[SZ_FNAME]
+int	beam_stat[ARB], len_flat
+
+int	i
+int	beam, len[MAX_NR_BEAMS]
+char	sfname[SZ_FNAME]
+pointer	flatsp[MAX_NR_BEAMS], imflat
+
+int	strlen(), imgeti()
+pointer	imgl1r(), immap()
+errchk	immap
+
+begin
+	# Determine beam number.
+
+	iferr (beam = imgeti (im, "BEAM-NUM"))
+	    beam = 0
+	beam = beam + 1
+
+	# Validate beam number
+	if (beam < 1 || beam > MAX_NR_BEAMS) {
+	    call eprintf (" Beam number out of range: %d - using 0\n")
+		call pargi (beam)
+	    beam = 1
+	}
+
+	# Has this beam already been loaded?
+	if (IS_INDEFI (beam_stat[beam])) {
+
+	    # Create file name
+	    call strcpy (flat_file, sfname, SZ_FNAME)
+
+	    # Flat field file names have beam number appended
+	    call sprintf (sfname[strlen(sfname)+1], SZ_FNAME, ".%04d")
+		call pargi (beam-1)
+
+	    # Open spectrum
+	    imflat = immap (sfname, READ_ONLY, 0)
+
+	    # Allocate space for this beam's sensitivity spectrum
+	    call salloc (flatsp[beam], IM_LEN(imflat,1), TY_REAL)
+
+	    # Copy pixels into space
+	    call amovr (Memr[imgl1r(imflat)], Memr[flatsp[beam]],
+		IM_LEN(imflat,1))
+
+	    # Must be careful that no division by zero occurs.
+	    do i = 1, IM_LEN(imflat,1) {
+		if (Memr[flatsp[beam]+i-1] == 0.0)
+		    Memr[flatsp[beam]+i-1] = 1.0
+	    }
+
+	    # Mark this beam accounted for
+	    beam_stat[beam] = 1
+	    len[beam] = IM_LEN(imflat,1)
+
+	    call imunmap (imflat)
+	}
+
+	# Point to the spectrum
+	sp = flatsp[beam]
+	len_flat = len[beam]
+
+end
+
+# DIVIDE -- Perform the division and create new spectrum
+
+procedure divide (im, flat, len_flat, ifile, ofile, rec, 
+	coincidence, ccmode, dtime, power)
+
+pointer	im, flat
+int	len_flat, rec, ccmode
+real	dtime, power
+char	ifile[ARB], ofile[ARB]
+bool	coincidence
+
+real	itm, imgetr()
+int	i, co_flag, imgeti()
+int	ncols, nlines
+char	calfname[SZ_FNAME], original[SZ_FNAME]
+pointer	imcal, rawpix, calpix
+
+pointer	immap(), impl2r(), imgl2r()
+
+begin
+	# Find smallest length of the two possible spectra
+	ncols = min (IM_LEN (im, 1), len_flat)
+
+	# Create new spectrum.  Make up a name
+	call sprintf (calfname, SZ_FNAME, "%s.%04d")
+	    call pargstr (ofile)
+	    call pargi (rec)
+
+	call xt_mkimtemp (ifile, calfname, original, SZ_FNAME)
+	imcal = immap (calfname, NEW_COPY, im)
+
+	IM_NDIM(imcal) = IM_NDIM(im)
+	IM_LEN (imcal,1) = ncols
+	IM_PIXTYPE(imcal) = TY_REAL
+
+	# Check for 2D spectrum
+	if (IM_NDIM(im) > 1)
+	    nlines = IM_LEN(im,2)
+	else
+	    nlines = 1
+
+	# Copy across the image title
+	call strcpy (IM_TITLE(im), IM_TITLE(imcal), SZ_LINE)
+
+	# Operate on the pixels
+	do i = 1, nlines {
+	    rawpix = imgl2r (im,i)
+	    calpix = impl2r (imcal,i)
+
+	    # Apply coincidence correction if needed
+	    co_flag = -1
+	    if (coincidence) {
+		iferr (co_flag = imgeti (im, "CO-FLAG"))
+		    ;
+	        if (co_flag < 1) {
+		    iferr (itm = imgetr (im, "EXPOSURE"))
+			iferr (itm = imgetr (im, "ITIME"))
+			    iferr (itm = imgetr (im, "EXPTIME"))
+				itm = 1.
+		    call coincor (Memr[rawpix], Memr[rawpix], ncols,
+			itm, co_flag, dtime, power, ccmode)
+		}
+	    }
+
+	    call adivr (Memr[rawpix], Memr[flat], Memr[calpix], ncols)
+	}
+
+	call imaddi (imcal, "QD-FLAG", 0)
+	if (co_flag != -1)
+	    call imaddi (imcal, "CO-FLAG", co_flag)
+
+	# Send user report
+	call printf ("writing [%s]: %s\n")
+	    call pargstr (original)
+	    call pargstr (IM_TITLE(imcal))
+	call flush (STDOUT)
+
+	call imunmap (im)
+	call imunmap (imcal)
+	call xt_delimtemp (calfname, original)
+end
+
+# INIT_BS -- Initialize beam status flags
+
+procedure init_bs (beam_stat)
+
+int	beam_stat[ARB]
+
+int	i
+
+begin
+	do i = 1, MAX_NR_BEAMS
+	    beam_stat[i] = INDEFI
+end
diff --git a/noao/onedspec/irsiids/t_flatfit.x b/noao/onedspec/irsiids/t_flatfit.x
new file mode 100644
index 00000000..c3558afc
--- /dev/null
+++ b/noao/onedspec/irsiids/t_flatfit.x
@@ -0,0 +1,740 @@
+include	<imhdr.h>
+include	<math/curfit.h>
+include	<gset.h>
+
+define	MAX_NR_BEAMS	100	# Max number of instrument apertures
+
+define	KEY		"noao$lib/scr/flatfit.key"
+define	PROMPT		"flatfit cursor options"
+
+# Definitions for Plotting modes
+define	PLT_FIT	1		# Plot the direct fit
+define	PLT_ERR	2		# Plot the errors in the fit
+define	PLT_LIN	3		# Plot the fit minus the linear part
+
+# T_FLATFIT -- Accumulate a series of flat field spectra to produce
+#  a grand sum and fit a function to the sum to produce a normalized
+#  flat containing the pixel-to-pixel variations.
+#  User interaction via the graphics cursor is provided. The following
+#  cursor commands are recognized:
+#
+#	? - Screen help
+#	/ - Status line help
+#	e - Plot in residual error mode
+#	f - Plot in fit to the data mode
+#	o - Change order of fit
+#	l - Change lower rejection sigma
+#	u - Change upper rejection sigma
+#	r - Reset fit to include rejected pixels
+#	s - Change upper and lower sigmas to same value
+#	i - Iterate again
+#	n - Iterate N times
+#	q - Quit and accept current solution (also RETURN)
+#
+
+procedure t_flatfit ()
+
+pointer	image					# Image name to be fit
+pointer	images					# Image name to be fit
+pointer	ofile					# Output image file name
+int	function				# Fitting function
+int	order					# Order of fitting function
+int	recs					# Spectral record numbers
+int	root, nrecs				# CL and ranges flags
+real	expo					# Exposure time
+real	dtime					# Deadtime
+real	power					# Power law coin. correction
+real	lower					# Lower rejection sigma
+real	upper					# Upper threshold sigma
+int	ngrow					# Rejection radius
+real	div_min					# Division min for option RESP
+bool	coincidence, all			# Apply coincidence correction
+bool	interact				# Interactive levels
+pointer	bstat					# Status of each aperture
+pointer	npts					# Length of spectrum
+pointer	esum					# Accumulated exposure time
+pointer	accum					# Pointers to beam accumulators
+pointer	title
+int	ccmode, beam
+int	niter
+
+int	i
+pointer	sp, str, im
+
+int	clgeti(), clgwrd(), clpopni(), imgeti()
+int	get_next_image(), decode_ranges()
+real	clgetr(), imgetr()
+bool	clgetb()
+pointer	immap()
+
+begin
+	call smark (sp)
+	call salloc (image, SZ_FNAME, TY_CHAR)
+	call salloc (images, MAX_NR_BEAMS, TY_POINTER)
+	call salloc (ofile, SZ_FNAME, TY_CHAR)
+	call salloc (recs, 300, TY_INT)
+	call salloc (bstat, MAX_NR_BEAMS, TY_INT)
+	call salloc (npts, MAX_NR_BEAMS, TY_INT)
+	call salloc (esum, MAX_NR_BEAMS, TY_REAL)
+	call salloc (accum, MAX_NR_BEAMS, TY_POINTER)
+	call salloc (title, MAX_NR_BEAMS, TY_POINTER)
+	call salloc (str, SZ_LINE, TY_CHAR)
+	call amovki (NULL, Memi[images], MAX_NR_BEAMS)
+
+	# Get task parameters.
+	root = clpopni ("input")
+
+	# Get input record numbers
+	call clgstr ("records", Memc[str], SZ_LINE)
+	if (decode_ranges (Memc[str], Memi[recs], 100, nrecs) == ERR)
+	    call error (0, "Bad range specification")
+
+	call clgstr ("output", Memc[ofile], SZ_LINE)
+
+	call clgcurfit ("function", "order", function, order)
+
+	lower = clgetr ("lower")
+	upper = clgetr ("upper")
+	ngrow = clgeti ("ngrow")
+	div_min = clgetr ("div_min")
+
+	# Determine desired level of activity
+	interact = clgetb ("interact")
+	all      = clgetb ("all_interact")
+
+	niter = clgeti ("niter")
+
+	# Is coincidence correction to be performed?
+	coincidence = clgetb ("coincor")
+
+	if (coincidence) {
+	    ccmode = clgwrd ("ccmode", Memc[str], SZ_LINE, ",photo,iids,")
+	    dtime  = clgetr ("deadtime")
+	    power = clgetr ("power")
+	}
+
+	call reset_next_image ()
+
+	# Clear all beam status flags
+	call amovki (INDEFI, Memi[bstat], MAX_NR_BEAMS)
+	call aclrr (Memr[esum], MAX_NR_BEAMS)
+
+	call printf ("Accumulating spectra --\n")
+	call flush (STDOUT)
+
+10	while (get_next_image (root, Memi[recs], nrecs, Memc[image],
+	    SZ_FNAME) != EOF) {
+	    iferr (im = immap (Memc[image], READ_ONLY, 0)) {
+		call eprintf ("Header info not available for [%s]\n")
+		    call pargstr (Memc[image])
+		goto 10
+	    }
+
+	    iferr (beam = imgeti (im, "BEAM-NUM"))
+		beam = 0
+	    if (beam < 0  || beam > MAX_NR_BEAMS-1)
+		call error (0, "Invalid aperture number")
+
+	    iferr (expo = imgetr (im, "EXPOSURE"))
+		iferr (expo = imgetr (im, "ITIME"))
+		    iferr (expo = imgetr (im, "EXPTIME"))
+			expo = 1
+
+	    # Add spectrum into accumulator
+	    if (IS_INDEFI (Memi[bstat+beam])) {
+	        Memi[npts+beam] = IM_LEN (im,1)
+	        call salloc (Memi[accum+beam], Memi[npts+beam], TY_REAL)
+	        call aclrr (Memr[Memi[accum+beam]], Memi[npts+beam])
+	        Memi[bstat+beam] = 0
+
+	        call salloc (Memi[title+beam], SZ_LINE, TY_CHAR)
+	        call strcpy (IM_TITLE(im), Memc[Memi[title+beam]], SZ_LINE)
+	    }
+
+	    call ff_accum_spec (im, Memi[npts], expo, Memi[bstat], beam+1,
+		Memi[accum], Memr[esum], coincidence, ccmode, dtime, power,
+		Memi[title])
+
+	    call printf ("[%s] added to aperture %1d\n")
+		call pargstr (Memc[image])
+		call pargi (beam)
+	    call flush (STDOUT)
+	    if (Memi[images+beam] == NULL)
+		call salloc (Memi[images+beam], SZ_FNAME, TY_CHAR)
+	    call strcpy (Memc[image], Memc[Memi[images+beam]], SZ_FNAME)
+
+	    call imunmap (im)
+	}
+
+	# Review all apertures containing data and perform fits.
+	# Act interactively if desired
+	do i = 0, MAX_NR_BEAMS-1 {
+	    if (!IS_INDEFI (Memi[bstat+i])) {
+		call fit_spec (Memr[Memi[accum+i]], Memi[npts+i], Memr[esum+i],
+		    interact, function, order, niter, lower, upper, ngrow,
+		    div_min, i)
+		if (interact & !all)
+		    interact = false
+		call wrt_fit_spec (Memc[Memi[images+i]], Memr[Memi[accum+i]],
+		    Memr[esum+i], Memc[ofile], i, Memc[Memi[title+i]],
+		    Memi[npts+i], order)
+	    }
+	}
+
+	call sfree (sp)
+	call clpcls (root)
+end
+
+# ACCUM_SPEC -- Accumulate spectra by beams
+
+procedure ff_accum_spec (im, len, expo, beam_stat, beam, accum, expo_sum,
+	coincidence, ccmode, dtime, power, title)
+
+pointer	im, accum[ARB], title[ARB]
+real	expo, expo_sum[ARB]
+int	beam_stat[ARB], beam, len[ARB]
+bool	coincidence
+int	ccmode
+real	dtime, power
+
+int	npts, co_flag, imgeti()
+pointer	pix
+
+pointer	imgl1r()
+
+begin
+	npts = IM_LEN (im, 1)
+
+	# Map pixels and optionally correct for coincidence
+	pix = imgl1r (im)
+	if (coincidence) {
+	    iferr (co_flag = imgeti (im, "CO-FLAG"))
+		co_flag = -1
+	    if (co_flag < 1) {
+		call coincor (Memr[pix], Memr[pix], npts, expo, co_flag,
+		    dtime, power, ccmode)
+	    }
+	}
+
+	# Add in the current data
+	npts = min (npts, len[beam])
+
+	call aaddr (Memr[pix], Memr[accum[beam]], Memr[accum[beam]], npts)
+
+	beam_stat[beam] = beam_stat[beam] + 1
+	expo_sum [beam] = expo_sum [beam] + expo
+end
+
+# WRT_FIT_SPEC -- Write out normalized spectrum
+
+procedure wrt_fit_spec (image, accum, expo_sum, ofile, beam, title, npts, order)
+
+char	image[SZ_FNAME]
+real	accum[ARB], expo_sum
+int	beam, npts, order
+char	ofile[SZ_FNAME]
+char	title[SZ_LINE]
+
+char	output[SZ_FNAME], temp[SZ_LINE]
+pointer	im, imnew, newpix
+
+pointer	immap(), impl1r()
+int	strlen()
+
+begin
+	im = immap (image, READ_ONLY, 0)
+10	call strcpy (ofile, output, SZ_FNAME)
+	call sprintf (output[strlen (output) + 1], SZ_FNAME, ".%04d")
+	    call pargi (beam)
+
+	# Create new image with a user area
+	# If an error occurs, ask user for another name to try
+	# since many open errors result from trying to overwrite an
+	# existing image.
+
+	iferr (imnew = immap (output, NEW_COPY, im)) {
+	    call eprintf ("Cannot create [%s] -- Already exists??\07\n")
+		call pargstr (output)
+	    call clgstr ("output", ofile, SZ_FNAME)
+	    go to 10
+	}
+
+	call strcpy ("Normalized flat:", temp, SZ_LINE)
+	call sprintf (temp[strlen (temp) + 1], SZ_LINE, "%s")
+	    call pargstr (title)
+	call strcpy (temp, IM_TITLE (imnew), SZ_LINE)
+	IM_PIXTYPE (imnew) = TY_REAL
+
+	newpix = impl1r (imnew)
+	call amovr (accum, Memr[newpix], npts)
+
+	call imaddr (imnew, "EXPOSURE", expo_sum)
+	call imaddi (imnew, "QF-FLAG", order)
+	call imunmap (im)
+	call imunmap (imnew)
+
+	call printf ("Fit for aperture %1d --> [%s]\n")
+	    call pargi (beam)
+	    call pargstr (output)
+	call flush (STDOUT)
+end
+
+# FIT_SPEC -- Fit a line through the spectrum with user interaction
+
+procedure fit_spec (accum, npts, expo_sum, interact, function, 
+	order, niter, lower, upper, ngrow, div_min, beam)
+
+real	accum[ARB], expo_sum
+bool	interact
+int	function, order, niter, ngrow, npts, beam
+real	lower, upper, div_min
+
+int	cc, key, gp, plt_mode
+int	i, initer, sum_niter, newgraph
+real	x1, y1, sigma, temp
+pointer	sp, wts, x, y, cv
+bool	first
+char	gtitle[SZ_LINE], command[SZ_FNAME]
+
+int	clgcur(), clgeti()
+pointer	gopen()
+real	clgetr(), cveval()
+
+data	plt_mode/PLT_FIT/
+
+begin
+	# Perform initial fit
+	call smark (sp)
+	call salloc (wts, npts, TY_REAL)
+	call salloc (x  , npts, TY_REAL)
+	call salloc (y  , npts, TY_REAL)
+
+	first = true
+	if (!interact) {
+	    sum_niter = 0
+	    do i = 1, niter
+		call linefit (accum, npts, function, order, lower, upper, 
+		    ngrow, cv, first, Memr[wts], Memr[x])
+	    sum_niter = niter
+
+	} else {
+	    gp = gopen ("stdgraph", NEW_FILE, STDGRAPH)
+	    call sprintf (gtitle, SZ_LINE, "Flat Field Sum - %f seconds ap:%1d")
+		call pargr (expo_sum)
+		call pargi (beam)
+
+	    key = 'r'
+	    repeat {
+		switch (key) {
+		    case 'e': # Plot errors
+			plt_mode = PLT_ERR
+			newgraph = YES
+
+		    case 'f': # Plot fit
+			plt_mode = PLT_FIT
+			newgraph = YES
+
+		    case 'o': # Change order
+			order = clgeti ("new_order")
+			# Reinstate all pixels
+			first = true
+			newgraph = YES
+
+		    case 'l': # Change lower sigma
+			lower = clgetr ("new_lower")
+			newgraph = YES
+
+		    case 'u': # Change upper sigma
+			upper = clgetr ("new_upper")
+			newgraph = YES
+
+		    case 'r': # Reset fit parameters
+			first = true
+			newgraph = YES
+
+		    case 's': # Change both rejection sigmas
+			lower = clgetr ("new_lower")
+			upper = lower
+			call clputr ("new_upper", upper)
+			newgraph = YES
+
+		    case 'i': # Iterate again - Drop thru
+			initer = 1
+			newgraph = YES
+
+		    case 'n': # Iterate n times
+			initer = clgeti ("new_niter")
+			newgraph = YES
+
+		    case 'q': # Quit
+			break
+
+		    case '?': # Clear and help
+			call gpagefile (gp, KEY, PROMPT)
+
+		    case '/': # Status line help
+			call ff_sts_help
+
+		    case 'I': # Interrupt
+			call fatal (0, "Interrupt")
+
+		    default:
+			call printf ("\07\n")
+		}
+
+		if (newgraph == YES) {
+		    # Suppress an iteration if plot mode change requested
+		    if (key != 'e' && key != 'f') {
+		        if (first) {
+			    sum_niter = 0
+			    initer = niter
+			    call cvfree (cv)
+		        }
+		        do i = 1, initer
+			    call linefit (accum, npts, function, order, lower, 
+			        upper, ngrow, cv, first, Memr[wts], Memr[x])
+		        sum_niter = sum_niter + initer
+		    }
+
+		    switch (plt_mode) {
+		    case PLT_FIT:
+		        call plot_fit (gp, accum, cv, function, order, npts, 
+			    gtitle, Memr[wts], Memr[x], Memr[y], sigma)
+		    case PLT_ERR:
+		        call plot_fit_er (gp, accum, cv, function, order, npts,
+			    gtitle, Memr[wts], Memr[x], Memr[y], sigma)
+		    }
+
+		    newgraph = NO
+	        }
+	    } until (clgcur ("cursor",x1,y1,cc,key,command,SZ_FNAME) == EOF)
+	    call gclose (gp)
+	}
+
+	# Replace original data with the data/fit
+	do i = 1, npts {
+	    temp = cveval (cv, real (i))
+	    if (temp == 0.0)
+		temp = max (temp, div_min)
+	    accum[i] = accum[i] / temp
+	}
+
+	call cvfree (cv)
+	call sfree (sp)
+
+	# Save iteration count for next time
+	niter = sum_niter
+end
+
+# LINEFIT -- Fit desired function thru data
+
+procedure linefit (pix, npts, function, order, lower, upper, ngrow, cv,
+	    first, wts, x)
+
+real	pix[ARB]			# Data array to fit
+int	npts				# Elements in array
+int	function			# Type of fitting function
+int	order				# Order of fitting function
+real	lower				# Lower rejection threshold
+real	upper				# Upper rejection threshold
+int	ngrow				# Rejection growing radius
+pointer	cv
+real	wts[ARB]			# Array weights
+real	x[ARB]
+bool	first
+
+int	ier, i, nreject
+
+int	reject()
+
+begin
+10	if (first) {
+	    do i = 1, npts {
+		x[i] = i
+		wts[i] = 1.0
+	    }
+
+	    # Initialize curve fitting.
+	    call cvinit (cv, function, order, 1., real (npts))
+	    call cvfit (cv, x, pix, wts, npts, WTS_USER, ier)
+	    nreject = 0
+	    first = false
+	}
+
+	# Do pixel rejection if desired.
+	if ((lower > 0.) || (upper > 0.))
+	        nreject = reject (cv, x, pix, wts, npts, lower, upper, ngrow)
+	else
+		nreject = 0
+
+	if (nreject == ERR) {
+	    call eprintf ("Cannot fit data -- too many points rejected??\n")
+	    call cvfree (cv)
+	    first = true
+	    go to 10
+	}
+end
+
+# REJECT -- Reject points with large residuals from the fit.
+#
+# The sigma of the input to the fit is calculated.  The rejection thresholds
+# are set at -lower*sigma and upper*sigma.  Points outside the rejection
+# thresholds are rejected from the fit and flaged by setting their
+# weights to zero.  Finally, the remaining points are refit and a new
+# fit line evaluated.  The number of points rejected is returned.
+
+int procedure reject (cv, x, y, w, npoints, lower, upper, ngrow)
+
+pointer	cv				# Curve descriptor
+real	x[ARB]				# Input ordinates
+real	y[ARB]				# Input data values
+real	w[ARB]				# Weights
+int	npoints				# Number of input points
+real	lower				# Lower rejection sigma
+real	upper				# Upper rejection sigma
+int	ngrow				# Rejection radius
+
+int	i, j, n, i_min, i_max, nreject
+real	sigma, residual, resid_min, resid_max
+
+real	cveval()
+
+begin
+	# Determine sigma of fit and set rejection limits.
+	sigma = 0.
+	n = 0
+	do i = 1, npoints {
+	    if (w[i] == 0.)
+		next
+	    sigma = sigma + (y[i] - cveval (cv, x[i])) ** 2
+	    n = n + 1
+	}
+
+	sigma = sqrt (sigma / (n - 1))
+	resid_min = -lower * sigma
+	resid_max = upper * sigma
+
+	# Reject the residuals exceeding the rejection limits.
+	nreject = 0
+	for (i = 1; i <= npoints; i = i + 1) {
+	    if (w[i] == 0.)
+		next
+	    residual = y[i] - cveval (cv, x[i])
+	    if ((residual < resid_min) || (residual > resid_max)) {
+		i_min = max (1, i - ngrow)
+		i_max = min (npoints, i + ngrow)
+
+		# Reject points from the fit and flag them with zero weight.
+		do j = i_min, i_max {
+		    call cvrject (cv, x[j], y[j], w[j])
+		    w[j] = 0.
+		    nreject = nreject + 1
+		}
+		i = i_max
+	    }
+	}
+
+	# Refit if points have been rejected.
+	if (nreject > 0) {
+	    call cvsolve (cv, i)
+	    if (i != OK)
+		return (ERR)
+	}
+
+	return (nreject)
+end
+
+# PLOT_FIT -- Plot the fit to the image line and data
+
+procedure plot_fit (gp, pix, cv, function, order, npts, gtitle, wts, xfit,
+	yfit, sigma)
+
+int	gp, npts, function, order
+real	pix[ARB], wts[ARB],  xfit[ARB], yfit[ARB]
+pointer	cv
+real	sigma
+char	gtitle[SZ_LINE]
+
+real	x1, x2
+int	i
+
+begin
+	# Set up plot
+	x1 = 1.0
+	x2 = npts
+
+	call gseti (gp, G_NMINOR, 0)
+	call gclear (gp)
+	call gsview (gp, 0.15, 0.95, 0.20, 0.9)
+	call gploto (gp, pix, npts, x1, x2, gtitle)
+
+	# Now plot the fit
+	do i = 1, npts
+	    xfit[i] = i
+
+	call cvvector (cv, xfit, yfit, npts)
+	call gvline (gp, yfit, npts, x1, x2)
+
+	# Compute sigma and write it out
+	call get_sigma (pix, yfit, wts, npts, sigma)
+	call show_status (function, order, sigma, npts, wts)
+end
+
+# PLOT_FIT_ER -- Plot the error in the fit to the image line and data
+
+procedure plot_fit_er (gp, pix, cv, function, order, npts, gtitle, wts, xfit,
+	yfit, sigma)
+
+int	gp, npts, function, order
+real	pix[ARB], wts[ARB], xfit[ARB], yfit[ARB]
+pointer	cv
+real	sigma
+char	gtitle[SZ_LINE]
+
+real	x1, x2, y[2]
+int	i
+
+begin
+	# Set up plot
+	x1 = 1.0
+	x2 = npts
+	y[1] = -0.0001
+	y[2] = +0.0001
+
+	call cvvector (cv, xfit, yfit, npts)
+
+	# Compute percentage errors
+	do i = 1, npts
+	    if (pix[i] != 0.0)
+		yfit[i] = (pix[i] - yfit[i]) / pix[i]
+	    else
+		yfit[i] = 0.0
+
+	call gseti (gp, G_NMINOR, 0)
+	call gclear (gp)
+	call gsview (gp, 0.15, 0.95, 0.20, 0.9)
+
+	call gploto (gp, yfit, npts, x1, x2, 
+	   "Flat field fractional error in fit")
+
+	# Draw a zero error line
+	call gline (gp, x1, y[1], x2, y[2])
+
+	# Compute sigma
+	call get_sigma0 (yfit, wts, npts, sigma)
+	call show_status (function, order, sigma, npts, wts)
+end
+
+# SHOW_STATUS -- Show the fit status on status line
+
+procedure show_status (function, order, sigma, npts, wts)
+
+int	function, order, npts
+real 	sigma, wts[ARB]
+
+int	i, nvals
+
+begin
+	# Count non-rejected points
+	nvals = 0
+	do i = 1, npts
+	    if (wts[i] != 0.0)
+		nvals = nvals + 1
+
+	call printf ("Fit type: %s    order: %2d    rms: %6.3f")
+	switch (function) {
+	    case LEGENDRE:
+		call pargstr ("Legendre")
+	    case CHEBYSHEV:
+		call pargstr ("Chebyshev")
+	    case SPLINE3:
+		call pargstr ("Spline3")
+	    case SPLINE1:
+		call pargstr ("Spline1")
+	    default:
+		call pargstr ("???")
+	}
+
+	    call pargi (order)
+	    call pargr (sigma)
+
+	call printf ("   points: %d out of %d")
+	    call pargi (nvals)
+	    call pargi (npts)
+
+	call flush (STDOUT)
+end
+
+# GET_SIGMA -- Compute rms error between two vectors whose average difference
+#              is zero.
+
+procedure get_sigma (y1, y2, wts, n, sigma)
+
+real	y1[ARB], y2[ARB], wts[ARB], sigma
+int	n
+
+int	i, nval
+real	sum
+
+begin
+	sum = 0.0
+	nval = 0
+	do i = 1, n
+	    if (wts[i] != 0.0) {
+		sum = sum + (y1[i] - y2[i]) ** 2
+		nval = nval + 1
+	    }
+
+	sigma = sqrt (sum / (nval-1))
+	return
+end
+
+# GET_SIGMA0 -- Compute rms error of a vector
+
+procedure get_sigma0 (y1, wts, n, sigma)
+
+real	y1[ARB], wts[ARB], sigma
+int	n
+
+int	i, nval
+real	sum
+
+begin
+	sum = 0.0
+	nval = 0
+	do i = 1, n
+	    if (wts[i] != 0.0) {
+		sum = sum + y1[i]**2
+		nval = nval + 1
+	    }
+
+	sigma = sqrt (sum / (nval-1))
+	return
+end
+
+# FF_STS_HELP -- Status line help for Flat Fit
+
+procedure ff_sts_help ()
+
+int	linenr, maxline
+
+data	linenr/1/
+data	maxline/2/
+
+begin
+	switch (linenr) {
+	    case 1:
+	    call printf ("e=err plot  f=data plot  o=order  l=lower sigma  ")
+	    call printf ("u=upper sigma  s=both sigmas")
+
+	    case 2:
+	    call printf ("r=incl reject  i=iterate  n=niterate  q=quit  ")
+	    call printf ("?=help  /=linehelp  <CR>=quit")
+	}
+
+	call flush (STDOUT)
+
+	linenr = linenr + 1
+	if (linenr > maxline)
+	    linenr = 1
+end
diff --git a/noao/onedspec/irsiids/t_slist1d.x b/noao/onedspec/irsiids/t_slist1d.x
new file mode 100644
index 00000000..75837d50
--- /dev/null
+++ b/noao/onedspec/irsiids/t_slist1d.x
@@ -0,0 +1,163 @@
+include	<error.h>
+include	<imhdr.h>
+include	<fset.h>
+include <smw.h>
+
+
+# SLIST1D --  Lists header information from IIDS/IRS format header
+# This is the original T_SLIST.
+
+procedure t_slist1d ()
+
+int	root
+int	long_header
+pointer	sp, image, im, mw, sh, ptr
+int	i, nl, df, sm, qf, qd, bs, co
+
+int	btoi(), imtgetim(), imgeti()
+bool	clgetb()
+pointer	imtopenp(), immap(), smw_openim()
+errchk	immap, smw_openim, shdr_open
+
+begin
+	call smark (sp)
+	call salloc (image, SZ_LINE, TY_CHAR)
+
+	# Parameters
+	root = imtopenp ("input")
+	call clgstr ("records", Memc[image], SZ_LINE)
+	call odr_openp (root, Memc[image])
+	long_header = btoi (clgetb ("long_header"))
+
+	# Initialize
+	call fseti (STDOUT, F_FLUSHNL, YES)
+
+	# Loop over all input images by subsets
+	while (imtgetim (root, Memc[image], SZ_FNAME) != EOF) {
+
+	    # Open image
+	    iferr {
+		im = NULL
+		mw = NULL
+		ptr = immap (Memc[image], READ_ONLY, 0); im = ptr
+		ptr = smw_openim (im); mw = ptr
+		call shdr_open (im, mw, 1, 1, INDEFI, SHHDR, sh)
+	    } then {
+		if (mw != NULL) {
+		    call smw_close (mw)
+		    if (sh != NULL)
+			MW(sh) = NULL
+		}
+		if (im != NULL)
+		    call imunmap (im)
+		call erract (EA_WARN)
+		next
+	    }
+
+	    nl = IM_LEN(im,2)
+	    do i = 1, nl {
+		call shdr_open (im, mw, i, 1, INDEFI, SHHDR, sh)
+
+		if (long_header == YES) {
+		    call printf ("[%s] %4dpts %s\n")
+			call pargstr (IMNAME(sh))
+			call pargi (SN(sh))
+			call pargstr (TITLE(sh))
+
+		    if (OFLAG(sh) == 1) {
+			call printf ("oflag = OBJECT, beam_number = %d")
+			    call pargi (BEAM(sh))
+		    } else if (OFLAG (sh) == 0) {
+			call printf ("oflag = SKY,    beam_number = %d")
+			    call pargi (BEAM(sh))
+		    }
+		    call printf (",\n")
+
+		    iferr (df = imgeti (im, "DF-FLAG"))
+			df = -1
+		    iferr (sm = imgeti (im, "SM-FLAG"))
+			sm = -1
+		    iferr (qf = imgeti (im, "QF-FLAG"))
+			qf = -1
+		    iferr (qd = imgeti (im, "QD-FLAG"))
+			qd = -1
+		    iferr (bs = imgeti (im, "BS-FLAG"))
+			bs = -1
+		    iferr (co = imgeti (im, "CO-FLAG"))
+			co = -1
+
+		    # Airmass may not be in header. It could be computed if
+		    # if the observatory latitude were available.
+
+		    call printf ("airmass = %5.3f,%25tW0 = %0.3f,")
+			call pargr (AM(sh))
+			call pargr (W0(sh))
+		    call printf ("   WPC = %0.5g,      ITM = %.2f,\n")
+			call pargr (WP(sh))
+			call pargr (IT(sh))
+		    call printf ("NP1 = %d, NP2 = %d,")
+			call pargi (NP1(sh))
+			call pargi (NP2(sh))
+		    call printf ("    UT = %0.1h,   ST = %0.1h,\n")
+			call pargr (UT(sh))
+			call pargr (ST(sh))
+		    call printf ("HA = %0.2h,")
+			call pargr (HA(sh))
+		    call printf ("        RA = %0.2h,   DEC = %0.1h,\n")
+			call pargr (RA(sh))
+			call pargr (DEC(sh))
+		    call printf (
+			"df = %d, sm = %d, qf = %d, dc = %d, qd = %d, ") 
+			call pargi (df)
+			call pargi (sm)
+			call pargi (qf)
+			call pargi (DC(sh))
+			call pargi (qd)
+		    call printf ("ex = %d, bs = %d, ca = %d, co = %d")
+			call pargi (EC(sh))
+			call pargi (bs)
+			call pargi (FC(sh))
+			call pargi (co)
+
+		    call printf ("\n\n")
+		} else {
+		    if (nl == 1) {
+			call printf ("[%s]:%s %4ds %4dpts   %s\n")
+			    call pargstr (IMNAME(sh))
+			    if (OFLAG(sh) == 1)
+				call pargstr ("o")
+			    else
+				call pargstr ("s")
+			    call pargr (IT(sh))
+			    call pargi (SN(sh))
+			    call pargstr (TITLE(sh))
+		    } else {
+			call printf ("[%s]:%s %6.2fs %4dpts %dspectra   %s\n")
+			    call pargstr (IMNAME(sh))
+			    if (OFLAG(sh) == 1)
+				call pargstr ("o")
+			    else
+				call pargstr ("s")
+			    call pargr (IT(sh))
+			    call pargi (SN(sh))
+			    call pargi (nl)
+			    call pargstr (TITLE(sh))
+		    }
+		    break
+		}
+	    }
+
+	    call smw_close (mw)
+	    if (sh != NULL)
+		MW(sh) = NULL
+	    call imunmap (im)
+	}
+
+	# Null out record string to avoid learn mode
+	call clpstr ("records", "")
+
+	# Free space
+	call shdr_close (sh)
+	call imtclose (root)
+	call sfree (sp)
+end
diff --git a/noao/onedspec/irsiids/t_subsets.x b/noao/onedspec/irsiids/t_subsets.x
new file mode 100644
index 00000000..6a2a61bf
--- /dev/null
+++ b/noao/onedspec/irsiids/t_subsets.x
@@ -0,0 +1,121 @@
+include	<error.h>
+include	<imhdr.h>
+
+
+# T_SUBSETS -- Sub a series of spectra by pairs. A single spectrum
+#  is produced for every pair.
+#
+
+procedure t_subsets ()
+
+pointer	image
+pointer	recstr, ofile
+int	root, start_rec, subset
+int	nrecs
+int	npts, nrem, ifile, tog
+real	expo, wtsum
+pointer	sp, recs, im[2], cur_pix, sp_sum
+
+real	imgetr()
+int	clpopni(), clgeti()
+int	get_next_image(), decode_ranges()
+pointer	immap(), imgl1r()
+
+begin
+	call smark (sp)
+	call salloc (image, SZ_FNAME, TY_CHAR)
+	call salloc (ofile, SZ_FNAME, TY_CHAR)
+	call salloc (recstr, SZ_LINE, TY_CHAR)
+	call salloc (recs, 300, TY_INT)
+
+	# Open input file name template
+	root   = clpopni ("input")
+
+	# Get range specification if any
+	call clgstr ("records", Memc[recstr], SZ_LINE)
+	if (decode_ranges (Memc[recstr], Memi[recs], 100, nrecs) == ERR)
+	    call error (0, "Bad range specification")
+
+	# Get rootname for output files and starting record
+	call clgstr ("output", Memc[ofile], SZ_FNAME)
+	start_rec = clgeti ("start_rec")
+
+	# Initialize range decoder
+	call reset_next_image ()
+
+	#Initialize file counter
+	ifile = 0
+
+	# Set weighting value needed by spectrum writer
+	wtsum = 1.0
+
+	# Define subset of operation is a pair
+	subset = 2
+
+	# Loop over all input images by subsets
+	while (get_next_image (root, Memi[recs], nrecs, Memc[image], 
+	    SZ_FNAME) != EOF) {
+
+	    # Get toggle value
+	    tog = mod (ifile, 2) + 1
+
+	    # Open image
+	    iferr (im[tog] = immap (Memc[image], READ_ONLY, 0)) {
+		call erract (EA_WARN)
+		next
+	    }
+
+	    # Load data
+	    cur_pix = imgl1r (im[tog])
+
+	    # Allocate space for the sum
+	    if (mod (ifile,2) == 0) {
+		npts = IM_LEN (im[tog],1)
+		call calloc (sp_sum, npts, TY_REAL)
+
+		# Zero exposure counter
+		expo = 0.0
+
+		# Add first spectrum
+	        call amovr (Memr[cur_pix], Memr[sp_sum], npts)
+
+		iferr (expo = imgetr (im[tog], "EXPOSURE"))
+		    iferr (expo = imgetr (im[tog], "ITIME"))
+			iferr (expo = imgetr (im[tog], "EXPTIME"))
+			    expo = 1
+
+	        call printf ("[%s] added\n")
+		    call pargstr (Memc[image])
+		call flush (STDOUT)
+
+	    } else {
+		# Subtract second spectrum
+		call asubr (Memr[sp_sum], Memr[cur_pix], Memr[sp_sum], 
+		    min (npts, IM_LEN(im[tog],1)))
+		call printf ("[%s] subtracted\n")
+		    call pargstr (Memc[image])
+		call flush (STDOUT)
+		call imunmap (im[2])
+
+		call wrt_set (Memr[sp_sum], subset, im[1], Memc[ofile],
+		    start_rec, expo, wtsum, -1)
+		call mfree (sp_sum, TY_REAL)
+	    }
+
+	    ifile = ifile + 1
+	}
+	# Check that there are no remaining spectra in an unfulfilled subset
+	nrem = mod (ifile, 2)
+	if (nrem != 0) {
+	    call mfree (sp_sum, TY_REAL)
+
+	    call eprintf ("Unfulfilled pair ignored\n")
+	}
+
+	# Update record number
+	call clputi ("next_rec", start_rec)
+
+	# Free space
+	call sfree (sp)
+	call clpcls (root)
+end
diff --git a/noao/onedspec/irsiids/t_sums.x b/noao/onedspec/irsiids/t_sums.x
new file mode 100644
index 00000000..e28ebb35
--- /dev/null
+++ b/noao/onedspec/irsiids/t_sums.x
@@ -0,0 +1,239 @@
+include	<error.h>
+include	<imhdr.h>
+
+define	MAX_NR_BEAMS	100	# Max number of instrument apertures
+
+# T_SUMS -- Compute sums of strings of spectra according to
+#	    Aperture number and object/sky flag. So for IIDS/IRS
+#           type spectra, 4 sums will be generated.
+#           In general, there will be 2N sums where N is the number
+#           apertures.
+
+procedure t_sums ()
+
+pointer	image					# Image name to be added
+pointer	images					# Image name to be added
+pointer	ofile					# Output image file name
+pointer	recstr					# Record number string
+int	recs					# Spectral record numbers
+int	root, nrecs				# CL and ranges flags
+real	expo					# Exposure time
+pointer	bstat[2]				# Status of each aperture
+pointer	npts[2]					# Length of spectrum
+pointer	esum[2]					# Accumulated exposure time
+pointer	accum[2]				# Pointers to beam accumulators
+pointer	title[2]
+int	beam, object
+int	start_rec
+
+int	i, j
+pointer	sp, work, im
+
+real	imgetr()
+int	clgeti(), clpopni(), imgeti()
+int	get_next_image(), decode_ranges()
+pointer	immap()
+
+begin
+	call smark (sp)
+	call salloc (image, SZ_FNAME, TY_CHAR)
+	call salloc (images, MAX_NR_BEAMS, TY_POINTER)
+	call salloc (ofile, SZ_FNAME, TY_CHAR)
+	call salloc (recstr, SZ_LINE, TY_CHAR)
+	call salloc (recs, 300, TY_INT)
+	call salloc (accum, MAX_NR_BEAMS, TY_POINTER)
+	call salloc (title, MAX_NR_BEAMS, TY_POINTER)
+	call amovki (NULL, Memi[images], MAX_NR_BEAMS)
+	call salloc (work, 2*5*MAX_NR_BEAMS, TY_STRUCT)
+	bstat[1] = work
+	bstat[2] = work + MAX_NR_BEAMS
+	npts[1] = work + 2 * MAX_NR_BEAMS
+	npts[2] = work + 3 * MAX_NR_BEAMS
+	esum[1] = work + 4 * MAX_NR_BEAMS
+	esum[2] = work + 5 * MAX_NR_BEAMS
+	accum[1] = work + 6 * MAX_NR_BEAMS
+	accum[2] = work + 7 * MAX_NR_BEAMS
+	title[1] = work + 8 * MAX_NR_BEAMS
+	title[2] = work + 9 * MAX_NR_BEAMS
+
+	# Get task parameters.
+	root = clpopni ("input")
+
+	# Get input record numbers
+	call clgstr ("records", Memc[recstr], SZ_LINE)
+	if (decode_ranges (Memc[recstr], Memi[recs], 100, nrecs) == ERR)
+	    call error (0, "Bad range specification")
+
+	call clgstr ("output", Memc[ofile], SZ_LINE)
+
+	start_rec = clgeti ("start_rec")
+
+	call reset_next_image ()
+
+	# Clear all beam status flags
+	call amovki (INDEFI, Memi[bstat[1]], MAX_NR_BEAMS*2)
+	call aclrr (Memr[esum[1]], MAX_NR_BEAMS*2)
+
+	call printf ("Accumulating spectra --\n")
+	call flush (STDOUT)
+
+	while (get_next_image (root, Memi[recs], nrecs, Memc[image],
+	    SZ_FNAME) != EOF) {
+	    iferr (im = immap (Memc[image], READ_ONLY, 0)) {
+		call erract (EA_WARN)
+		next
+	    }
+
+	    # Load header
+	    iferr (beam = imgeti (im, "BEAM-NUM"))
+		beam = 0
+	    if (beam < 0  || beam > MAX_NR_BEAMS-1)
+		call error (0, "Invalid aperture number")
+
+	    # Select array: Object = array 2; sky = array 1
+	    iferr (object = imgeti (im, "OFLAG"))
+		object = 1
+	    if (object == 1)
+		object = 2
+	    else
+		object = 1
+
+	    iferr (expo = imgetr (im, "EXPOSURE"))
+		iferr (expo = imgetr (im, "ITIME"))
+		    iferr (expo = imgetr (im, "EXPTIME"))
+			expo = 1
+
+	    # Add spectrum into accumulator
+	    if (IS_INDEFI (Memi[bstat[object]+beam])) {
+	        Memi[npts[object]+beam] = IM_LEN (im,1)
+	        call salloc (Memi[accum[object]+beam], IM_LEN(im,1), TY_REAL)
+	        call aclrr (Memr[Memi[accum[object]+beam]], IM_LEN(im,1))
+	        Memi[bstat[object]+beam] = 0
+
+	        call salloc (Memi[title[object]+beam], SZ_LINE, TY_CHAR)
+	        call strcpy (IM_TITLE(im), Memc[Memi[title[object]+beam]],
+		    SZ_LINE)
+	    }
+
+	    call su_accum_spec (im, Memi[npts[1]], expo, Memi[bstat[1]],
+		beam+1, Memi[accum[1]], Memr[esum[1]], Memi[title[1]], object)
+
+	    call printf ("[%s] %s spectrum added to aperture %1d\n")
+		call pargstr (Memc[image])
+		if (object == 2)
+		    call pargstr ("object")
+		else
+		    call pargstr ("sky   ")
+		call pargi (beam)
+	    call flush (STDOUT)
+
+	    if (Memi[images+beam] == NULL)
+		call salloc (Memi[images+beam], SZ_FNAME, TY_CHAR)
+	    call strcpy (Memc[image], Memc[Memi[images+beam]], SZ_FNAME)
+	    call imunmap (im)
+	}
+
+	# Review all apertures containing data and write sums
+	do i = 0, MAX_NR_BEAMS-1
+	    do j = 1, 2
+	    if (!IS_INDEFI (Memi[bstat[j]+i])) {
+		call wrt_spec (Memc[Memi[images+i]], Memr[Memi[accum[j]+i]],
+		    Memr[esum[j]+i], Memc[ofile], start_rec,
+		    Memc[Memi[title[j]+i]], Memi[npts[j]+i], i, j)
+
+		start_rec = start_rec + 1
+	    }
+
+	call clputi ("next_rec", start_rec)
+	call sfree (sp)
+	call clpcls (root)
+end
+
+# ACCUM_SPEC -- Accumulate spectra by beams
+
+procedure su_accum_spec (im, len, expo, beam_stat, beam, accum, expo_sum,
+	title, object)
+
+pointer	im, accum[MAX_NR_BEAMS,2], title[MAX_NR_BEAMS,2]
+real	expo, expo_sum[MAX_NR_BEAMS,2]
+int	beam_stat[MAX_NR_BEAMS,2], beam, len[MAX_NR_BEAMS,2]
+int	object
+
+int	npts
+pointer	pix
+
+pointer	imgl1r()
+
+begin
+	npts = IM_LEN (im, 1)
+
+	# Map pixels and optionally correct for coincidence
+	pix = imgl1r (im)
+
+	# Add in the current data
+	npts = min (npts, len[beam, object])
+
+	call aaddr (Memr[pix], Memr[accum[beam, object]], 
+	    Memr[accum[beam, object]], npts)
+
+	beam_stat[beam, object] = beam_stat[beam, object] + 1
+	expo_sum [beam, object] = expo_sum [beam, object] + expo
+end
+
+# WRT_SPEC -- Write out normalized spectrum
+
+procedure wrt_spec (image, accum, expo_sum, ofile, start, title, npts, object,
+	beam)
+
+char	image[SZ_FNAME]
+real	accum[ARB], expo_sum
+int	start, npts
+char	ofile[SZ_FNAME]
+char	title[SZ_LINE]
+int	object, beam
+
+char	output[SZ_FNAME], temp[SZ_LINE]
+pointer	im, imnew, newpix
+
+pointer	immap(), impl1r()
+int	strlen()
+
+begin
+	im = immap (image, READ_ONLY, 0)
+10	call strcpy (ofile, output, SZ_FNAME)
+	call sprintf (output[strlen (output) + 1], SZ_FNAME, ".%04d")
+	    call pargi (start)
+
+	# Create new image with a user area
+	# If an error occurs, ask user for another name to try
+	# since many open errors result from trying to overwrite an
+	# existing image.
+
+	iferr (imnew = immap (output, NEW_COPY, im)) {
+	    call eprintf ("Cannot create [%s] -- Already exists??\07\n")
+		call pargstr (output)
+	    call clgstr ("newoutput", ofile, SZ_FNAME)
+	    go to 10
+	}
+
+	call strcpy ("Summation:", temp, SZ_LINE)
+	call sprintf (temp[strlen (temp) + 1], SZ_LINE, "%s")
+	    call pargstr (title)
+	call strcpy (temp, IM_TITLE (imnew), SZ_LINE)
+
+	newpix = impl1r (imnew)
+	call amovr (accum, Memr[newpix], npts)
+
+	call imaddr (imnew, "EXPOSURE", expo_sum)
+	call imunmap (im)
+	call imunmap (imnew)
+
+	call printf ("%s sum for aperture %1d --> [%s]\n")
+	    if (object == 1)
+		call pargstr ("Object")
+	    else
+		call pargstr ("Sky   ")
+	    call pargi (beam)
+	    call pargstr (output)
+	call flush (STDOUT)
+end
diff --git a/noao/onedspec/irsiids/t_widstape.x b/noao/onedspec/irsiids/t_widstape.x
new file mode 100644
index 00000000..1f96d146
--- /dev/null
+++ b/noao/onedspec/irsiids/t_widstape.x
@@ -0,0 +1,343 @@
+include	<mach.h>
+include	<error.h>
+include	<imhdr.h>
+include <smw.h>
+
+define	SZ_IDSTITLE	64		# Length of IDSOUT title
+define	SZ_CARD		80		# Columns on a card
+
+# T_WIDSTAPE -- Convert each line of an IRAF image to IDSOUT text format.
+# Each image line is treated as a one dimensional spectrum.
+# A maximum IDSOUT length of 1024 points is enforced silently.
+#
+# There are two types of output:
+#    single   -- All image lines are appended to a single IDSOUT file.
+#    multiple -- Each image line is appended to a different IDSOUT file.
+
+procedure t_widstape ()
+
+pointer	image			# Image to be converted
+pointer	recs			# Record numbers
+pointer	idsout			# IDSOUT file or root name to be written
+int	block_size		# Block size
+bool	ebcdic			# ASCII or EBCDIC
+
+int	i, mfd, root, nrecs
+pointer	sp, im, mw, sh, ptr
+
+int	open(), mtopen(), clgeti(), clpopni()
+int	get_next_image(), decode_ranges(), mtfile(), mtneedfileno()
+bool	clgetb()
+pointer	immap(), smw_openim()
+errchk	immap, smw_openim, shdr_open, wrt_ids_rec
+
+begin
+	call smark (sp)
+	call salloc (image, SZ_LINE, TY_CHAR)
+	call salloc (idsout, SZ_FNAME, TY_CHAR)
+	call salloc (recs, 300, TY_INT)
+
+	# Parameters
+	root   = clpopni ("input")
+	call clgstr ("records", Memc[image], SZ_LINE)
+	call clgstr ("idsout", Memc[idsout], SZ_FNAME)
+	block_size = clgeti ("block_size")
+	ebcdic = clgetb ("ebcdic")
+
+	# Set record numbers
+	if (decode_ranges (Memc[image], Memi[recs], 100, nrecs) == ERR)
+	    call error (0, "Bad range specification")
+
+	# Check that a realistic block size was requested
+	if (mod (block_size, SZ_CARD) == 0)
+	    block_size = block_size / SZB_CHAR
+	else
+	    call error (0, "Blocks not integral number of cards")
+
+	# Open output tape file
+	# First determine if a file number was specified
+	if (mtfile (Memc[idsout]) == YES) {
+
+	    # If no file, check if new_tape was specified and if so,
+	    # force file=1; otherwise force file=EOT
+
+	    if (mtneedfileno (Memc[idsout]) == YES) {
+		if (!clgetb("new_tape"))
+		    call mtfname (Memc[idsout], EOT, Memc[idsout], SZ_FNAME)
+
+		else
+		    call mtfname (Memc[idsout], 1, Memc[idsout], SZ_FNAME)
+	    }
+	    mfd = mtopen (Memc[idsout], WRITE_ONLY, block_size)
+	} else
+	    mfd = open (Memc[idsout], NEW_FILE, BINARY_FILE)
+
+	# Loop over all files
+	call reset_next_image ()
+	while (get_next_image (root, Memi[recs], nrecs, Memc[image],
+	    SZ_LINE) != EOF) {
+	    iferr {
+		im = NULL
+		mw = NULL
+		ptr = immap (Memc[image], READ_ONLY, 0); im = ptr
+		ptr = smw_openim (im); mw = ptr
+
+		# Write out a spectrum for each line in the image
+		do i = 1, IM_LEN (im,2) {
+		    call shdr_open (im, mw, i, 1, INDEFI, SHDATA, sh)
+		    call wrt_ids_rec (mfd, sh, Memc[image], ebcdic)
+		}
+
+		call printf ("copied - [%s]: %s\n")
+		call pargstr (IMNAME(sh))
+		call pargstr (TITLE(sh))
+		call flush (STDOUT)
+	    } then
+		call erract (EA_WARN)
+
+	    if (mw != NULL)
+		call smw_close (mw)
+	    if (im != NULL)
+		call imunmap (im)
+	}
+
+	call shdr_close (sh)
+	call close (mfd)
+	call sfree (sp)
+end
+
+
+# WRT_IDS_REC -- Write one IIDS/IRS format record in IDSOUT form
+
+procedure wrt_ids_rec (mfd, sh, image, ebcdic)
+
+int	mfd
+pointer	sh
+char	image[SZ_FNAME]
+bool	ebcdic
+
+# IDSOUT header parameters
+char	label[SZ_IDSTITLE]		# Record label
+int	record				# Record number
+int	uttime				# UT time in seconds
+int	st				# Siderial time in seconds
+real	ra				# Right Ascension in hours
+real	dec				# Declination in degrees
+real	ha				# Hour angle in hours
+real	airmass				# Air mass
+int	itime				# Integration time
+real	wavelen1			# Wavelength of first pixel
+real	dispersion			# Dispersion per pixel
+
+int     i, rec_no, df, sm, qf, qd, bs, co
+pointer	sp, padline, bufline
+
+int	strmatch(), imgeti()
+
+begin
+	call smark (sp)
+	call salloc (padline, SZ_LINE, TY_CHAR)
+	call salloc (bufline, SZ_LINE, TY_CHAR)
+
+	# Fill in header parameters.
+
+	call strcpy (TITLE(sh), label, SZ_IDSTITLE)
+
+	# The following two calculations were causing floating overflows
+	# when the header values were indefinite.  SEH 7-23-86
+	if (IS_INDEF(UT(sh)))
+	    uttime = INDEFI
+	else
+	    uttime = UT(sh) * 3600.
+
+	if (IS_INDEF(ST(sh)))
+	    st = INDEFI
+	else
+	    st     = ST(sh) * 3600.
+
+	ra         = RA(sh)
+	dec        = DEC(sh)
+	ha         = HA(sh)
+	airmass    = AM(sh)
+	itime      = IT(sh)
+	wavelen1   = W0(sh)
+	dispersion = WP(sh)
+
+	iferr (df = imgeti (IM(sh), "DF-FLAG"))
+	    df = -1
+	iferr (sm = imgeti (IM(sh), "SM-FLAG"))
+	    sm = -1
+	iferr (qf = imgeti (IM(sh), "QF-FLAG"))
+	    qf = -1
+	iferr (qd = imgeti (IM(sh), "QD-FLAG"))
+	    qd = -1
+	iferr (bs = imgeti (IM(sh), "BS-FLAG"))
+	    bs = -1
+	iferr (co = imgeti (IM(sh), "CO-FLAG"))
+	    co = -1
+
+	# Create a padding line to fill the IDSOUT block to 1024 points.
+
+	call sprintf (Memc[padline], SZ_LINE,
+	    "%10.4e%10.4e%10.4e%10.4e%10.4e%10.4e%10.4e%10.4e\n")
+	do i = 1, 8
+	    call pargr (0.)
+
+	# Line 1 -- Record number, etc.
+	rec_no = strmatch (image, ".")
+	call sscan (image[rec_no])
+	    call gargi (record)
+
+	    call sprintf (Memc[bufline], SZ_LINE, 
+		"%5d%5d%15.7e%15.7e%5d%5d%5d%5d%5d%5d%10d")
+		call pargi (record)
+		call pargi (itime)
+		call pargr (wavelen1)
+		call pargr (dispersion)
+		call pargi (0)
+		call pargi (SN(sh))
+		call pargi (BEAM(sh))
+		call pargi (-1)
+		call pargi (-1)
+		call pargi (0)
+		call pargi (uttime)
+
+	    call putcard (mfd, Memc[bufline], ebcdic)
+
+	    # Line 2 -- Siderial time, RA, and Dec.
+
+	    call sprintf (Memc[bufline], SZ_LINE, 
+		"%10d%15.7e%15.7e%5d%5d%5d%5d%5d%5d%5d%5d")
+		call pargi (st)
+		call pargr (ra)
+		call pargr (dec)
+		call pargi (0)
+		call pargi (df)
+		call pargi (sm)
+		call pargi (qf)
+		call pargi (DC(sh))
+		call pargi (qd)
+		call pargi (EC(sh))
+		call pargi (bs)
+
+	    call putcard (mfd, Memc[bufline], ebcdic)
+
+	    # Line 3 -- Hour angle, air mass, UT date, and exposure title.
+
+	    call sprintf (Memc[bufline], SZ_LINE, 
+		"%5d%5d%2w%-3.3s%5d%15.7e%15.7e%27wEND")
+		call pargi (FC(sh))
+		call pargi (co)
+		call pargstr ("IRF")
+		call pargi (OFLAG(sh))
+		call pargr (ha)
+		call pargr (airmass)
+
+	    call putcard (mfd, Memc[bufline], ebcdic)
+
+	    # Line 4 -- Record label.
+	    call sprintf (Memc[bufline], SZ_LINE, "%-77sEND")
+		call pargstr (TITLE(sh))
+
+	    call putcard (mfd, Memc[bufline], ebcdic)
+
+	    # Lines 5 to 132
+
+	    call putdata (mfd, Memr[SY(sh)], SN(sh), Memc[padline],
+		Memc[bufline], ebcdic)
+
+	    # Line 133 -- Blank line
+
+	    call sprintf (Memc[bufline], SZ_LINE, "%80w")
+	    call putcard (mfd, Memc[bufline], ebcdic)
+end
+
+
+# PUTDATA -- Format and output extraction data to IDSOUT length of 1024 points.
+# Special effort is made to make the zero padding efficient.
+
+procedure putdata (mfd, data, npts, padline, bufline, ebcdic)
+
+int	mfd				# IDSOUT file descriptor
+real	data[npts]			# Data
+int	npts				# Number of data points
+char	padline[ARB]			# Padding string
+char	bufline[ARB]			# Output buffer string
+bool	ebcdic				# Convert to ebcdic
+
+int	i, j, k, l, n
+int	index
+double	ddata
+
+int	dtoc3()
+
+begin
+	j = min (1024, npts)	# Maximum number of data points
+	k = j / 8 * 8		# Index of last data point in last complete line
+	if (k < j)
+	   l = k + 8		# Index of last point in last line with data
+	else
+	   l = k
+
+	# Write all complete data lines.
+
+	index = 1
+	do i = 1, k {
+	    ddata = double (data[i])
+	    n = dtoc3 (ddata, bufline[index], 10, 4, 'e', 10)
+	    while (n < 10) {
+		bufline[index+n] = ' '
+		n = n + 1
+	    }
+	    index = index + 10
+	    if (mod (i, 8) == 0) {
+	        call putcard (mfd, bufline, ebcdic)
+		index = 1
+	    }
+	}
+
+	# Write partial data line.
+
+	index = 1
+	do i = k + 1, l {
+	    if (i <= j) {
+		ddata = double (data[i])
+	        n = dtoc3 (ddata, bufline[index], 11, 5, 'e', 10)
+	    } else
+	        n = dtoc3 (0.D0, bufline[index], 11, 5, 'e', 10)
+	    while (n < 10) {
+		bufline[index+n] = ' '
+		n = n + 1
+	    }
+	    index = index + 10
+	    if (mod (i, 8) == 0) {
+	        call putcard (mfd, bufline, ebcdic)
+		index = 1
+	    }
+	}
+
+	# Write remaining padding lines.
+
+	do i = l + 1, 1024, 8
+	    call putcard (mfd, padline, ebcdic)
+end
+
+# PUTCARD -- Convert to ebcdic if desired and write out card
+
+procedure putcard (mfd, bufline, ebcdic)
+
+int	mfd
+char	bufline[ARB]
+bool	ebcdic
+
+char	packline[SZ_LINE]
+
+begin
+	if (ebcdic) {
+	    call ascii_to_ebcdic (bufline, packline, SZ_CARD)
+	    call achtsb (packline, packline, SZ_CARD)
+	} else
+	    call chrpak (bufline, 1, packline, 1, SZ_CARD)
+
+	call write (mfd, packline, SZ_CARD/SZB_CHAR)
+end
diff --git a/noao/onedspec/irsiids/widstape.par b/noao/onedspec/irsiids/widstape.par
new file mode 100644
index 00000000..33dee906
--- /dev/null
+++ b/noao/onedspec/irsiids/widstape.par
@@ -0,0 +1,8 @@
+# IDSOUT parameter file -- write a CYBER style IDSOUT tape
+
+idsout,s,a,,,,Output file or magtape
+input,s,a,,,,Image root name to write
+records,s,a,,,,Records to write
+block_size,i,h,3200,80,10640,Tape block size in bytes
+new_tape,b,h,no,,,Is this a new (blank) tape
+ebcdic,b,h,no,,,Convert character code to ebcdic