From fa080de7afc95aa1c19a6e6fc0e0708ced2eadc4 Mon Sep 17 00:00:00 2001
From: Joseph Hunkeler <jhunkeler@gmail.com>
Date: Wed, 8 Jul 2015 20:46:52 -0400
Subject: Initial commit

---
 noao/onedspec/sensfunc/sfoutput.x | 114 ++++++++++++++++++++++++++++++++++++++
 1 file changed, 114 insertions(+)
 create mode 100644 noao/onedspec/sensfunc/sfoutput.x

(limited to 'noao/onedspec/sensfunc/sfoutput.x')

diff --git a/noao/onedspec/sensfunc/sfoutput.x b/noao/onedspec/sensfunc/sfoutput.x
new file mode 100644
index 00000000..e21df280
--- /dev/null
+++ b/noao/onedspec/sensfunc/sfoutput.x
@@ -0,0 +1,114 @@
+include	<mach.h>
+include	<imhdr.h>
+include	"sensfunc.h"
+
+
+# SF_OUTPUT -- Write the sensitivity function image.
+
+procedure sf_output (stds, nstds, cv, output, ignoreaps)
+
+pointer	stds[nstds]		# Standard star data
+int	nstds			# Number of standard stars
+pointer	cv			# Sensitivity function curve
+char	output[SZ_FNAME]	# Output root image name (must be SZ_FNAME)
+bool	ignoreaps		# Ignore apertures?
+
+int	i, ap, nw, scan(), nscan()
+real	w1, w2, dw, dw1, aplow[2], aphigh[2], cveval()
+pointer	sp, fname, std, im, mw, buf, immap(), mw_open(), impl1r()
+errchk	imaddi, imaddr
+
+define	makeim_	99
+
+begin
+	# Return if no output root sensitivity imagename is specified.
+	if (output[1] == EOS)
+	    return
+
+	call smark (sp)
+	call salloc (fname, SZ_FNAME, TY_CHAR)
+
+	# Determine wavelength range and reference standard star.
+	w1 = MAX_REAL
+	w2 = -MAX_REAL
+	dw = MAX_REAL
+	do i = 1, nstds-2 {
+	    if (STD_FLAG(stds[i]) == SF_EXCLUDE)
+		next
+	    std = stds[i]
+	    dw1 = abs ((STD_WEND(std) - STD_WSTART(std)) / (STD_NPTS(std) - 1))
+	    w1 = min (w1, STD_WSTART(std), STD_WEND(std))
+	    w2 = max (w2, STD_WSTART(std), STD_WEND(std))
+	    dw = min (dw, dw1)
+	}
+	nw = (w2 - w1) / dw + 1.5
+
+	# Make output image name with aperture number appended.  If the
+	# image exists allow the user to change root name.
+makeim_
+	if (ignoreaps) {
+	    call strcpy (output, Memc[fname], SZ_FNAME)
+	} else {
+	    call sprintf (Memc[fname], SZ_FNAME, "%s.%04d")
+		call pargstr (output)
+		call pargi (STD_BEAM(std))
+	}
+
+	iferr (im = immap (Memc[fname], NEW_IMAGE, 0)) {
+	    call printf ("Cannot create %s --  Enter new name: ")
+		call pargstr (Memc[fname])
+	    call flush (STDOUT)
+	    if (scan() != EOF) {
+	        call gargwrd (Memc[fname], SZ_FNAME)
+	        if (nscan() == 1) {
+	    	    call strcpy (Memc[fname], output, SZ_FNAME)
+	    	    goto makeim_
+		}
+	    }
+	    call printf ("No sensitivity function created for aperture %2d\n")
+	        call pargi (STD_BEAM(std))
+	    call flush (STDOUT)
+	    return
+	}
+
+	# Define the image header.
+	IM_NDIM(im) = 1
+	IM_LEN(im,1) = nw
+	IM_PIXTYPE(im) = TY_REAL
+	if (ignoreaps) {
+	    call sprintf (IM_TITLE(im), SZ_FNAME,
+		"Sensitivity function for all apertures")
+	} else {
+	    call sprintf (IM_TITLE(im), SZ_FNAME,
+		"Sensitivity function for aperture %d")
+		call pargi (STD_BEAM(std))
+	}
+
+	mw = mw_open (NULL, 1)
+	call mw_newsystem (mw, "equispec", 1)
+	call mw_swtype (mw, 1, 1, "linear", "label=Wavelength units=Angstroms")
+	call smw_open (mw, NULL, im)
+	ap = STD_BEAM(std)
+	aplow[1] = INDEF
+	aphigh[1] = INDEF
+	aplow[2] = INDEF
+	aphigh[2] = INDEF
+	call smw_swattrs (mw, 1, 1, ap, STD_BEAM(std), 0,
+	    double(w1), double(dw), nw, 0D0, aplow, aphigh, "")
+	call smw_saveim (mw, im)
+	call smw_close (mw)
+
+	# Write sensitivity data.
+	buf = impl1r (im)
+	do i = 0, nw-1
+	    Memr[buf+i] = cveval (cv, w1 + i * dw)
+
+	# Notify user.
+	call printf ("%s --> %s\n")
+		call pargstr (IM_TITLE(im))
+		call pargstr (Memc[fname])
+	call flush (STDOUT)
+
+	call imunmap (im)
+	call sfree (sp)
+end
-- 
cgit