path: root/noao/imred/src/fibers/proc.cl

# PROC -- Process spectra from 2D to wavelength calibrated 1D.
# This program combines the operations of extraction, flat fielding,
# fiber throughput correction, dispersion correction, and sky subtraction
# in as simple and noninteractive way as possible.  Certain assumptions
# are made about the data and the output.  A blank sky image, called a
# sky flat, may be used to determine the instrument throughput.  The data
# must all share the same position on the 2D image and the same
# dispersion solution apart from small instrumental changes which can be
# tracked automatically.

procedure proc (objects, apref, flat, throughput, arcs1, arcs2, arcreplace,
	arctable, fibers, apidtable, crval, cdelt, objaps, skyaps, arcaps,
	objbeams, skybeams, arcbeams, scattered, fitflat, recenter, edit,
	trace, arcap, clean, dispcor, savearcs, skyalign, skysubtract,
	skyedit, saveskys, splot, redo, update, batch, listonly)

string	objects			{prompt="List of object spectra"}

file	apref			{prompt="Aperture reference spectrum"}
file	flat			{prompt="Flat field spectrum"}
file	throughput		{prompt="Throughput file or image (optional)"}
string	arcs1			{prompt="List of arc spectra"}
string	arcs2			{prompt="List of shift arc spectra"}
file	arcreplace		{prompt="Special aperture replacements"}
file	arctable		{prompt="Arc assignment table (optional)\n"}

int	fibers			{prompt="Number of fibers"}
file	apidtable		{prompt="Aperture identifications"}
string	crval = "INDEF"		{prompt="Approximate wavelength"}
string	cdelt = "INDEF"		{prompt="Approximate dispersion"}
string	objaps			{prompt="Object apertures"}
string	skyaps			{prompt="Sky apertures"}
string	arcaps			{prompt="Arc apertures"}
string	objbeams		{prompt="Object beam numbers"}
string	skybeams		{prompt="Sky beam numbers"}
string	arcbeams		{prompt="Arc beam numbers\n"}

bool	scattered		{prompt="Subtract scattered light?"}
bool	fitflat			{prompt="Fit and ratio flat field spectrum?"}
bool	recenter		{prompt="Recenter object apertures?"}
bool	edit			{prompt="Edit/review object apertures?"}
bool	trace			{prompt="Trace object spectra?"}
bool	arcap			{prompt="Use object apertures for arcs?"}
bool	clean			{prompt="Detect and replace bad pixels?"}
bool	dispcor			{prompt="Dispersion correct spectra?"}
bool	savearcs		{prompt="Save internal arcs?"}
bool	skyalign		{prompt="Align sky lines?"}
bool	skysubtract		{prompt="Subtract sky?"}
bool	skyedit			{prompt="Edit the sky spectra?"}
bool	saveskys		{prompt="Save sky spectra?"}
bool	splot			{prompt="Plot the final spectrum?"}
bool	redo			{prompt="Redo operations if previously done?"}
bool	update			{prompt="Update spectra if cal data changes?"}
bool	batch			{prompt="Extract objects in batch?"}
bool	listonly		{prompt="List steps but don't process?\n"}

real	datamax = INDEF		{prompt="Max data value / cosmic ray threshold"}

string	ansskyedit = "yes"	{prompt="Edit the sky spectra?", mode="q"}
bool	newaps, newresp, newdisp, newarcs, dobatch

string	anssplot = "yes"	{prompt="Splot spectrum?", mode="q",
				 enum="no|yes|NO|YES"}

string	extn = ".ms"		{prompt="Extraction extension"}
struct	*fd1, *fd2, *fd3

begin
	string	imtype, mstype
	string	arcref1, arcref2, spec, arc, align=""
	string	arcref1ms, arcref2ms, specms, arcms, response
	string	objs, temp, temp1, done
	string	str1, str2, str3, str4, arcrefs, log1, log2
	bool	scat, reextract, extract, disp, disperr, sky, log
	bool	skyedit1, skyedit2, splot1, splot2
	int	i, j, n, nspec
	struct	err

	# Call a separate task to do the listing to minimize the size of
	# this script and improve it's readability.

	dobatch = no
	if (listonly) {
	    listonly (objects, apidtable, apref, flat, throughput, arcs1, arcs2,
		scattered, dispcor, skysubtract, redo, update)
	    bye
	}

	imtype = "." // envget ("imtype")
	i = stridx (",", imtype)
	if (i > 0)
	    imtype = substr (imtype, 1, i-1)
	mstype = ".ms" // imtype
	n = strlen (imtype)

	# Temporary files used repeatedly in this script.  Under some
	# abort circumstances these files may be left behind.

	objs = mktemp ("tmp$iraf")
	temp = mktemp ("tmp$iraf")
	temp1 = mktemp ("tmp$iraf")
	done = mktemp ("tmp$iraf")

	if (apidtable != "") {
	    i = strlen (apidtable)
	    if (i > n && substr (apidtable, i-n+1, i) == imtype) {
		apidtable = substr (apidtable, 1, i-n)
		i = strlen (apidtable)
	    }
	    for (j=1; j<=i && substr(apidtable,j,j)==" "; j+=1);
	    apidtable = substr (apidtable, j, i)
	}
	i = strlen (apidtable)
	if (i == 0)
	    extn = ".ms"
	else {
	    extn = apidtable
	    while (yes) {
		i = stridx ("/$]", extn)
		if (i == 0)
		    break
		j = strlen (extn)
		extn = substr (extn, i+1, j)
	    }
	    extn = extn // ".ms"
	}

	# Get query parameter.
	getspec (objects, objs)
	if (arctable == "" || arctable == " ") {
	    if (arcs2 == "" || arcs2 == " ")
	        arcrefs = arcs1
	    else
		arcrefs = arcs2
	} else
	    arcrefs = arctable
	arcref1 = ""
	arcref2 = ""

	# Rather than always have switches on the logfile and verbose flags
	# we use TEE and set a file to "dev$null" if output is not desired.
	# We must check for the null string to signify no logfile.

	tee.append = yes
	if (logfile == "")
	    log1 = "dev$null"
	else
	    log1 = logfile
	if (verbose)
	    log2 = "STDOUT"
	else
	    log2 = "dev$null"

	# If the update switch is used changes in the calibration data
	# can cause images to be reprocessed (if they are in the object
	# list).  Possible changes are in the aperture definitions,
	# response function, dispersion solution, and sensitivity
	# function.  The newarcs flag is used to only go through the arc
	# image headers once setting the reference spectrum, airmass, and
	# UT.

	newaps = no
	newresp = no
	newdisp = no
	newarcs = yes

	# Check if there are aperture definitions in the database and
	# define them if needed.  This is usually somewhat interactive.
	# Delete the database entry to start fresh if we enter this
	# because of a redo.  Set the newaps flag in case an update is
	# desired.

	i = strlen (apref)
	if (i > n && substr (apref, i-n+1, i) == imtype)
	    apref = substr (apref, 1, i-n)

	getspec (apref, temp)
	fd1 = temp
	if (fscan (fd1, apref) == EOF)
	    error (1, "No aperture reference")
	fd1 = ""; delete (temp, verify=no)

	# Initialize
	apscript.saturation = INDEF
	apscript.references = apref
	apscript.profiles = ""
	apscript.apidtable = apidtable
	apscript.nfind = fibers
	apscript.clean = clean
	if (splot) {
	    splot1 = yes
	    splot2 = yes
	} else {
	    splot1 = no
	    splot2 = no
	}

	reextract = redo
	if (reextract || !access (database // "/ap" // apref // extn)) {
	    if (!access (apref // imtype)) {
		printf ("Aperture reference spectrum not found - %s%s\n",
		    apref, imtype) | scan (err)
		error (1, err // "\nCheck settting of imtype")
	    }
	    print ("Set reference apertures for ", apref) | tee (log1)
	    if (access (database // "/ap" // apref))
		delete (database // "/ap" // apref, verify=no)
	    if (access (database // "/ap" // apref//extn))
		delete (database // "/ap" // apref // extn, verify=no)
	    apscript.ansresize = "yes"
	    apscript.ansedit = "yes"
	    apscript.ansfittrace = "yes"
	    apscript (apref, references="", ansfind="YES", ansrecenter="NO",
		anstrace="YES", ansextract="NO")
	    newaps = yes
	    copy (database//"/ap"//apref, database//"/ap"//apref//extn,
		verbose=no)
	} else {
	    if (access (database // "/ap" // apref))
		delete (database // "/ap" // apref, verify=no)
	    copy (database//"/ap"//apref//extn, database//"/ap"//apref,
		verbose=no)
	}

	if (recenter)
	    apscript.ansrecenter = "YES"
	else
	    apscript.ansrecenter = "NO"
	apscript.ansresize = "NO"
	if (edit)
	    apscript.ansedit = "yes"
	else
	    apscript.ansedit = "NO"
	if (trace)
	    apscript.anstrace = "YES"
	else
	    apscript.anstrace = "NO"
	if (scattered) {
	    apscript.ansfitscatter = "yes"
	    apscript.ansfitsmooth = "yes"
	}
	apscript.ansfittrace = "NO"
	apscript.ansextract = "YES"
	apscript.ansreview = "NO"
	if (skyedit) {
	    skyedit1 = yes
	    skyedit2 = yes
	} else {
	    skyedit1 = no
	    skyedit2 = no
	}

	# The next step is to setup the scattered light correction if needed.
	# We use the flat field image for the interactive setting unless
	# one is not used and then we use the aperture reference.
	# If these images have been  scattered light corrected we assume the
	# scattered light functions parameters are correctly set.

	i = strlen (flat)
	if (i > n && substr (flat, i-n+1, i) == imtype)
	    flat = substr (flat, 1, i-n)

	if (flat != "")
	    spec = flat
	else
	    spec = apref

	getspec (spec, temp)
	fd1 = temp
	if (fscan (fd1, spec) == EOF)
	    error (1, "No flat field")
	fd1 = ""; delete (temp, verify=no)

	scat = no
	if (scattered) {
	    if (redo && access (spec//"noscat"//imtype)) {
		imdelete (spec, verify=no)
		imrename (spec//"noscat", spec)
	    }
	    hselect (spec, "apscatte", yes) | scan (str1)
	    if (nscan() == 0)
		scat = yes
	}
	if (scat) {
	    print ("Subtract scattered light from ", spec) | tee (log1)
	    #apscript.ansfitscatter = "yes"
	    #apscript.ansfitsmooth = "yes"
	    imrename (spec, spec//"noscat")
	    apscript (spec//"noscat", output=spec, ansextract="NO",
		ansscat="YES", anssmooth="YES")
	    #apscript.ansfitscatter = "NO"
	    #apscript.ansfitsmooth = "NO"
	}

	# The next step is to process the flat field image which is used
	# as a flat field and a throughput correction.

	spec = throughput
	i = strlen (spec)
	if (i > n && substr (spec, i-n+1, i) == imtype)
	    spec = substr (spec, 1, i-n)
	specms = spec // mstype

	if (spec != "" && access (spec//imtype)) {
	    getspec (spec, temp)
	    fd1 = temp
	    if (fscan (fd1, spec) == EOF)
		error (1, "No flat field")
	    fd1 = ""; delete (temp, verify=no)

	    scat = no
	    if (scattered) {
		if (redo && access (spec//"noscat"//imtype)) {
		    imdelete (spec, verify=no)
		    imrename (spec//"noscat", spec)
		}
		hselect (spec, "apscatte", yes) | scan (str1)
		if (nscan() == 0)
		    scat = yes
	    }
	    if (scat) {
		print ("Subtract scattered light from ", spec) | tee (log1)
		imrename (spec, spec//"noscat")
		apscript (spec//"noscat", output=spec, ansextract="NO",
		    ansscat="YES", anssmooth="YES")
	    }
	}

	response = ""
	if (flat != "" || spec != "") {
	    if (extn == ".ms")
		response = flat // spec // "norm.ms"
	    else
		response = flat // spec // extn
	    reextract = redo || (update && newaps)
	    if (reextract || !access (response // imtype) || (redo && scat)) {
	        print ("Create response function ", response) | tee (log1)

	        if (access (response // imtype))
		    imdelete (response, verify=no)
	        if (access (flat //mstype))
		    imdelete (flat//mstype, verify=no)
	        if (access (specms))
		    imdelete (specms, verify=no)

	        fibresponse (flat, spec, apref, response, recenter=recenter,
		    edit=edit, trace=trace, clean=clean,
		    fitflat=fitflat, interactive=params.f_interactive,
		    function=params.f_function, order=params.f_order)

	        newresp = yes
	    }
	}

	# If not dispersion correcting we can go directly to extracting
	# the object spectra.  The reference arcs are the first on
	# the arc lists.  The processing of the reference arcs is done
	# by the task ARCREFS.

	if (dispcor) {
	    getspec (arcs1, temp)
	    fd1 = temp
	    if (fscan (fd1, arcref1) == EOF)
		error (1, "No reference arcs")
	    fd1 = ""; delete (temp, verify=no)
	    if (!access (arcref1 // imtype)) {
		printf ("Arc reference spectrum not found - %s%s\n",
		    arcref1, imtype) | scan (err)
		error (1, err // "\nCheck settting of imtype")
	    }
	    arcref1ms = arcref1 // extn
	    reextract = redo || (update && newaps)
	    if (reextract && access (arcref1ms//imtype))
	        imdelete (arcref1ms, verify=no)

	    getspec (arcs2, temp)
	    fd1 = temp
	    if (fscan (fd1, arcref2) == EOF)
		arcref2 = ""
	    else {
	        if (!access (arcref2 // imtype)) {
		    printf ("Arc reference spectrum not found - %s%s\n",
			arcref2, imtype) | scan (err)
		    error (1, err // "\nCheck settting of imtype")
		}
	        arcref2ms = arcref2 // extn
	        if (reextract && access (arcref2ms//imtype))
	            imdelete (arcref2ms, verify=no)
	    }
	    fd1 = ""; delete (temp, verify=no)

	    arcrefs (arcref1, arcref2, extn, arcreplace, apidtable, response,
		crval, cdelt, done, log1, log2)

	    # Define alignment if needed.
	    if (skyalign) {
		align = "align" // extn // imtype
		if (reextract)
		    imdelete (align, verify=no, >& "dev$null")
	    }
	}

	# Now we are ready to process the object spectra.

	reextract = redo || (update && (newaps || newresp || newdisp))
	fd1 = objs
	while (fscan (fd1, spec) != EOF) {
	    # Check if previously done; i.e. arc.
	    if (access (done)) {
	        fd2 = done
	        while (fscan (fd2, specms) != EOF)
		    if (spec == specms)
		        break
	        if (spec == specms)
		    next
	        fd2 = ""
	    }
	    if (!access (spec // imtype)) {
		printf ("Object spectrum not found - %s%s\n",
		    spec, imtype) | tee (log1)
		print ("Check settting of imtype")
	    }
	    specms = spec // mstype

	    # Determine required operations from the flags and image header.
	    scat = no
	    extract = no
	    disp = no
	    sky = no
	    if (scattered) {
		if (redo && access (spec//"noscat"//imtype)) {
		    imdelete (spec, verify=no)
		    imrename (spec//"noscat", spec)
		}
		hselect (spec, "apscatte", yes) | scan (str1)
		if (nscan() == 0)
		    scat = yes
	    }
	    if (reextract || !access (specms) || (redo && scat))
		extract = yes
	    else {
		hselect (specms, "dclog1", yes) | scan (str1)
		if (nscan () == 1) {
		    extract = update && newdisp
		    if (update && !newdisp)
		        # We really should check if REFSPEC will assign
			# different reference spectra.
			;
		} else
		    disp = dispcor

		hselect (specms, "skysub", yes) | scan (str1)
		if (nscan() == 0)
		    sky = skysubtract
	    }
		
	    if (extract) {
		disp = dispcor
		sky = skysubtract
	    }

	    # If fully processed go to the next object.
	    if (!extract && !disp && !sky)
		next

	    # If not interactive and the batch flag is set submit rest to batch.
	    if (batch && !skyedit1 && !skyedit2 && !splot1 && !splot2 &&
		apscript.ansedit == "NO" && (!skyalign || access (align))) {
		fd1 = ""; delete (objs, verify=no)
		apscript.ansfitscatter = "NO"
		apscript.ansfitsmooth = "NO"

		flprcache
		batch.objects = objects
		batch.datamax = datamax
		batch.response = response
		batch.arcs1 = arcs1
		batch.arcs2 = arcs2
		batch.arcref1 = arcref1
		batch.arcref2 = arcref2
		batch.arcreplace = arcreplace
		batch.apidtable = apidtable
		batch.arcrefs = arcrefs
		batch.extn = extn
		batch.objaps = objaps
		batch.skyaps = skyaps
		batch.arcaps = arcaps
		batch.objbeams = objbeams
		batch.skybeams = skybeams
		batch.arcbeams = arcbeams
		batch.done = done
		batch.logfile = log1
		batch.redo = reextract
		batch.update = update
		batch.scattered = scattered
		batch.arcap = arcap
		batch.dispcor = dispcor
		batch.savearcs = savearcs
		batch.skyalign = skyalign
		batch.skysubtract = skysubtract
		batch.saveskys = saveskys
		batch.newaps = newaps
		batch.newresp = newresp
		batch.newdisp = newdisp
		batch.newarcs = newarcs
		dobatch = yes
		return
	    }

	    # Process the spectrum in foreground.
	    if (extract) {
		if (access (specms))
		    imdelete (specms, verify=no)

		if (scat) {
		    print ("Subtract scattered light from ", spec) | tee (log1)
		    imrename (spec, spec//"noscat")
		    apscript (spec//"noscat", output=spec, ansextract="NO",
			ansscat="YES", anssmooth="YES")
		}

		print ("Extract object spectrum ", spec) | tee (log1)
		hselect (spec, "date-obs,ut,exptime", yes, > temp1)
		hselect (spec, "ra,dec,epoch,st", yes, >> temp1)
		fd3 = temp1
		if (fscan (fd3, str1, str2, str3) == 3) {
		    setjd (spec, observatory=observatory, date="date-obs",
			time="ut", exposure="exptime", jd="jd", hjd="",
			ljd="ljd", utdate=yes, uttime=yes, listonly=no,
			>> log1)
		    if (fscan (fd3, str1, str2, str3, str4) == 4)
			setairmass (spec, intype="beginning",
			    outtype="effective", exposure="exptime",
			    observatory=observatory, show=no, update=yes,
			    override=yes, >> log1)
		}
		fd3 = ""; delete (temp1, verify=no)
		apscript (spec, nsubaps=params.nsubaps, saturation=datamax)
		sapertures (specms, apertures="", apidtable=apidtable,
		    wcsreset=no, verbose=no, beam=INDEF, dtype=INDEF, w1=INDEF,
		    dw=INDEF, z=INDEF, aplow=INDEF, aphigh=INDEF, title=INDEF)
		if (response != "") {
		    if (params.nsubaps == 1)
			sarith (specms, "/", response, specms, w1=INDEF,
			    w2=INDEF, apertures="", bands="", beams="",
			    apmodulus=0, reverse=no, ignoreaps=no,
			    format="multispec", renumber=no, offset=0,
			    clobber=yes, merge=no, errval=0, verbose=no)
		    else {
			blkrep (response, temp, 1, params.nsubaps)
			sarith (specms, "/", temp, specms, w1=INDEF,
			    w2=INDEF, apertures="", bands="", beams="",
			    apmodulus=0, reverse=no, ignoreaps=yes,
			    format="multispec", renumber=no, offset=0,
			    clobber=yes, merge=no, errval=0, verbose=no)
			imdelete (temp, verify=no)
		    }
		}
	    }

	    disperr = no
	    if (disp) {
		# Fix arc headers if necessary.
		if (newarcs) {
		    getspec (arcs1, temp)
	    	    fd2 = temp
	    	    while (fscan (fd2, arc) != EOF) {
			hselect (arc, "date-obs,ut,exptime", yes, > temp1)
			hselect (arc, "ra,dec,epoch,st", yes, >> temp1)
			fd3 = temp1
			if (fscan (fd3, str1, str2, str3) == 3) {
			    setjd (arc, observatory=observatory, date="date-obs",
				time="ut", exposure="exptime", jd="jd", hjd="",
				ljd="ljd", utdate=yes, uttime=yes, listonly=no,
				>> log1)
			    if (fscan (fd3, str1, str2, str3, str4) == 4)
				setairmass (arc, intype="beginning",
				    outtype="effective", exposure="exptime",
				    observatory=observatory, show=no, update=yes,
				    override=yes, >> log1)
			}
			fd3 = ""; delete (temp1, verify=no)
	        	hedit (arc, "refspec1", arc, add=yes, verify=no,
		    	    show=no, update=yes)
	        	hedit (arc, "arctype", "henear", add=yes, verify=no,
		    	    show=no, update=yes)
	    	    }
	    	    fd2 = ""; delete (temp, verify=no)
		    getspec (arcs2, temp)
	    	    fd2 = temp
	    	    while (fscan (fd2, arc) != EOF) {
			hselect (arc, "date-obs,ut,exptime", yes, > temp1)
			hselect (arc, "ra,dec,epoch,st", yes, >> temp1)
			fd3 = temp1
			if (fscan (fd3, str1, str2, str3) == 3) {
			    setjd (arc, observatory=observatory,
				date="date-obs", time="ut", exposure="exptime",
				jd="jd", hjd="", ljd="ljd", utdate=yes,
				uttime=yes, listonly=no, >> log1)
			    if (fscan (fd3, str1, str2, str3, str4) == 4)
				setairmass (arc, intype="beginning",
				    outtype="effective", exposure="exptime",
				    observatory=observatory, show=no,
				    update=yes, override=yes, >> log1)

			}
			fd3 = ""; delete (temp1, verify=no)
	        	hedit (arc, "refspec1", arc, add=yes, verify=no,
		    	    show=no, update=yes)
	        	hedit (arc, "arctype", "shift", add=yes, verify=no,
		    	    show=no, update=yes)
	    	    }
	    	    fd2 = ""; delete (temp, verify=no)
		    newarcs = no
		}

		print ("Assign arc spectra for ", spec) | tee (log1)
		refspectra (spec, references=arcrefs,
		    apertures="", refaps="", ignoreaps=no,
		    select=params.select, sort=params.sort,
		    group=params.group, time=params.time,
		    timewrap=params.timewrap, override=yes, confirm=no,
		    assign=yes, logfiles="STDOUT", verbose=no) |
		    tee (log1, > log2)

		doarcs (spec, response, arcref1, arcref2, extn, arcreplace,
		    apidtable, arcaps, arcbeams, savearcs, reextract, arcap,
		    log1, no, done)

	        hselect (specms, "refspec1", yes, > temp)
	        fd2 = temp
	        i = fscan (fd2, arc)
	        fd2 = ""; delete (temp, verify=no)
		if (i < 1) {
		    print ("No arc reference assigned for ", spec) | tee (log1)
		    disperr = yes
		} else {
		    if (skyalign)
			doalign (spec, specms, align, arcref1ms, log1, no)
	            print ("Dispersion correct ", spec) | tee (log1)
		    dispcor (specms, "", linearize=params.linearize,
			database=database, table=arcref1ms, w1=INDEF,
			w2=INDEF, dw=INDEF, nw=INDEF, log=params.log,
			flux=params.flux, samedisp=no, global=no,
			ignoreaps=no, confirm=no, listonly=no,
			verbose=verbose, logfile=logfile)
		    if (params.nsubaps > 1) {
			imrename (specms, temp, verbose=no)
			scopy (temp, specms, w1=INDEF, w2=INDEF,
			    apertures="1-999", bands="", beams="", apmodulus=0,
			    offset=0, format="multispec", clobber=no, merge=no,
			    renumber=no, verbose=no)
			blkavg (temp, temp, 1, params.nsubaps, option="sum")
			imcopy (temp, specms//"[*,*]", verbose=no)
			imdelete (temp, verify=no)
		    }
		}
	    }

	    if (sky && !disperr) {
		str1 = ""
		if (skyaps != "")
		    str1 = "skyaps=" // skyaps
		if (skybeams != "")
		    str1 = str1 // " skybeams=" // skybeams
	        print ("Sky subtract ", spec, ": ", str1) | tee (log1)
		if (skyedit1) {
		    str1 = substr (ansskyedit, 1, 1)
		    if (str1 == "N" || str1 == "Y")
			skyedit1 = no
		    if (str1 == "n" || str1 == "N")
			skyedit2 = no
		    else
			skyedit2 = yes
		}
		skysub.reject = params.reject
		skysub (specms, output="", objaps=objaps, skyaps=skyaps,
		    objbeams=objbeams, skybeams=skybeams, skyedit=skyedit2,
		    combine=params.combine, scale=params.scale,
		    saveskys=saveskys, logfile=logfile)
		params.reject = skysub.reject
	        hedit (specms, "skysub", yes, add=yes, show=no, verify=no,
		    update=yes)
	    }

	    if (!disperr && (extract || disp || sky)) {
		if (splot1) {
		    print (specms, ":")
		    str1 = anssplot
		    if (str1 == "NO" || str1 == "YES")
			splot1 = no
		    if (str1 == "no" || str1 == "NO")
			splot2 = no
		    else
			splot2 = yes
		}
		if (splot2)
		    splot (specms)
	    }

	    print (spec, >> done)
	}
	fd1 = ""; delete (objs, verify=no)

	if (access (done))
	    delete (done, verify=no)
end