Repatch (from linux) of OSX IRAF

1 files changed, 136 insertions, 0 deletions

diff --git a/noao/digiphot/apphot/fitpsf/apsffit.x b/noao/digiphot/apphot/fitpsf/apsffit.x
new file mode 100644
index 00000000..2e6d6dbc
--- /dev/null
+++ b/noao/digiphot/apphot/fitpsf/apsffit.x
@@ -0,0 +1,136 @@
+include <mach.h>
+include "../lib/apphotdef.h"
+include "../lib/apphot.h"
+include "../lib/fitpsfdef.h"
+include "../lib/noisedef.h"
+include "../lib/fitpsf.h"
+
+# APSFFIT -- Procedure to fit an analytic function to the PSF.
+
+int procedure apsffit (ap, im, wx, wy)
+
+pointer	ap		# pointer to the apphot structure
+pointer	im		# pointer to the IRAF image
+real	wx, wy		# object coordinates
+
+int	ier, fier
+pointer	psf, nse
+real	datamin, datamax, dmin, dmax, threshold
+int	apfbuf(), apsfradgauss(), apsfelgauss(), apsfmoments()
+
+begin
+	# Initialize.
+	psf = AP_PPSF(ap)
+	nse = AP_NOISE(ap)
+	AP_PFXCUR(psf) = wx
+	AP_PFYCUR(psf) = wy
+	if (IS_INDEFR(wx) || IS_INDEFR(wy)) {
+            AP_OPFXCUR(psf) = INDEFR
+            AP_OPFYCUR(psf) = INDEFR
+	} else {
+            switch (AP_WCSOUT(ap)) {
+            case WCS_WORLD, WCS_PHYSICAL:
+                call ap_ltoo (ap, wx, wy, AP_OPFXCUR(psf), AP_OPFYCUR(psf), 1)
+            case WCS_TV:
+                call ap_ltov (im, wx, wy, AP_OPFXCUR(psf), AP_OPFYCUR(psf), 1)
+            default:
+                AP_OPFXCUR(psf) = wx
+                AP_OPFYCUR(psf) = wy
+            }
+	}
+	call amovkr (INDEFR, Memr[AP_PPARS(psf)], AP_MAXNPARS(psf))
+	call amovkr (INDEFR, Memr[AP_PPERRS(psf)], AP_MAXNPARS(psf))
+
+	# Fetch the buffer of pixels.
+	ier = apfbuf (ap, im, wx, wy)
+	if (ier == AP_NOPSFAREA)
+	    return (AP_NOPSFAREA)
+
+	# Compute the min and max of the data subraster.
+	if (IS_INDEFR(AP_DATAMIN(ap)))
+	    datamin = -MAX_REAL
+	else
+	    datamin = AP_DATAMIN(ap)
+	if (IS_INDEFR(AP_DATAMAX(ap)))
+	    datamax = MAX_REAL
+	else
+	    datamax = AP_DATAMAX(ap)
+
+	switch (AP_PSFUNCTION(psf)) {
+
+	case AP_RADGAUSS:
+
+	    fier = apsfradgauss (Memr[AP_PSFPIX(psf)], AP_PNX(psf), AP_PNY(psf),
+	        AP_POSITIVE(ap), AP_FWHMPSF(ap) * AP_SCALE(ap), datamin,
+		datamax, AP_NOISEFUNCTION(nse), AP_EPADU(nse),
+		AP_READNOISE(nse) / AP_EPADU(nse), AP_PMAXITER(psf),
+		AP_PK2(psf), AP_PNREJECT(psf), Memr[AP_PPARS(psf)],
+		Memr[AP_PPERRS(psf)], AP_PSFNPARS(psf))
+
+	    Memr[AP_PPARS(psf)+1] = Memr[AP_PPARS(psf)+1] + wx - AP_PXC(psf)
+	    Memr[AP_PPARS(psf)+2] = Memr[AP_PPARS(psf)+2] + wy - AP_PYC(psf) 
+
+	case AP_ELLGAUSS:
+
+	    fier = apsfelgauss (Memr[AP_PSFPIX(psf)], AP_PNX(psf), AP_PNY(psf),
+	        AP_POSITIVE(ap), AP_FWHMPSF(ap) * AP_SCALE(ap), datamin,
+		datamax, AP_NOISEFUNCTION(nse), AP_EPADU(nse),
+		AP_READNOISE(nse) / AP_EPADU(nse), AP_PMAXITER(psf),
+		AP_PK2(psf), AP_PNREJECT(psf), Memr[AP_PPARS(psf)],
+		Memr[AP_PPERRS(psf)], AP_PSFNPARS(psf))
+
+	    Memr[AP_PPARS(psf)+1] = Memr[AP_PPARS(psf)+1] + wx - AP_PXC(psf)
+	    Memr[AP_PPARS(psf)+2] = Memr[AP_PPARS(psf)+2] + wy - AP_PYC(psf) 
+
+	case AP_MOMENTS:
+
+	    call alimr (Memr[AP_PSFPIX(psf)], AP_PNX(psf) * AP_PNY(psf),
+	        dmin, dmax)
+	    dmin = max (dmin, datamin)
+	    dmax = min (dmax, datamax)
+	    threshold = 0.0
+
+	    if (AP_POSITIVE(ap) == YES)
+	        fier = apsfmoments (Memr[AP_PSFPIX(psf)], AP_PNX(psf),
+		    AP_PNY(psf), dmin + threshold, dmax,
+		    AP_POSITIVE(ap), Memr[AP_PPARS(psf)], Memr[AP_PPERRS(psf)],
+		    AP_PSFNPARS(psf))
+	    else
+	        fier = apsfmoments (Memr[AP_PSFPIX(psf)], AP_PNX(psf),
+		    AP_PNY(psf), dmax - threshold, dmin,
+		    AP_POSITIVE(ap), Memr[AP_PPARS(psf)],
+		    Memr[AP_PPERRS(psf)], AP_PSFNPARS(psf))
+
+	    Memr[AP_PPARS(psf)+1] = Memr[AP_PPARS(psf)+1] + wx - AP_PXC(psf)
+	    Memr[AP_PPARS(psf)+2] = Memr[AP_PPARS(psf)+2] + wy - AP_PYC(psf) 
+
+	default:
+
+	    # do nothing gracefully
+
+        }
+
+        switch (AP_WCSOUT(ap)) {
+        case WCS_WORLD, WCS_PHYSICAL:
+            call ap_ltoo (ap, Memr[AP_PPARS(psf)+1], Memr[AP_PPARS(psf)+2],
+		Memr[AP_PPARS(psf)+1], Memr[AP_PPARS(psf)+2], 1)
+        case WCS_TV:
+            call ap_ltov (im, Memr[AP_PPARS(psf)+1], Memr[AP_PPARS(psf)+2],
+		Memr[AP_PPARS(psf)+1], Memr[AP_PPARS(psf)+2], 1)
+        default:
+	    ;
+        }
+
+	# Return the appropriate error code.
+	if (fier == AP_OK) {
+	    if (ier == AP_PSF_OUTOFBOUNDS)
+		return (AP_PSF_OUTOFBOUNDS)
+	    else
+		return (AP_OK)
+	} else if (fier == AP_NPSF_TOO_SMALL) {
+	    call amovkr (INDEFR, Memr[AP_PPARS(psf)], AP_PSFNPARS(psf))
+	    call amovkr (INDEFR, Memr[AP_PPERRS(psf)], AP_PSFNPARS(psf))
+	    return (AP_NPSF_TOO_SMALL)
+	} else
+	    return (fier)
+end