Repatch (from linux) of OSX IRAF

1 files changed, 383 insertions, 0 deletions

diff --git a/noao/astutil/pdm/pdmcursor.x b/noao/astutil/pdm/pdmcursor.x
new file mode 100644
index 00000000..067ce966
--- /dev/null
+++ b/noao/astutil/pdm/pdmcursor.x
@@ -0,0 +1,383 @@
+include <mach.h>
+include <ctype.h>
+include <error.h>
+include	<gset.h>
+include	<fset.h>
+include	<math/curfit.h>
+include "pdm.h"
+
+define	PROMPT		"pdm options"
+
+# PDM_CURSOR -- Get the next command from the user in a graphics cursor loop.
+# Perform the requested function.
+
+procedure pdm_cursor (pdmp, ptype, infile, flip)
+
+pointer	pdmp			# pointer to PDM data structure
+int	ptype			# type of plot on the screen
+char	infile[SZ_FNAME]	# input file name
+bool	flip			# flip the y-axis scale
+
+real	xx, yy
+double	x, y			# cursor coordinates
+double	xmax, xmin
+double	period
+int	wcs			# wcs to which coordinates belong
+int	key			# keystroke value of cursor event
+int	stridxs(), gqverify(), ier
+int	clgcur(), i
+int	pdm_delete(), index
+int	pdm_undelete(), pdm_findmin()
+double	rx1, rx2
+double	dcveval()
+double	pdm_signif(), signif
+pointer	weights			# pointer to temporary weights array for icfit
+pointer	command			# string value, if any
+pointer	sp, sptemp
+errchk	pdm_colon, icg_fit, pdm_fitphase, pdm_alltheta
+errchk	pdm_phase, pdm_signif, pdm_ampep, pdm_findmin
+
+begin
+	call smark (sp)
+	call salloc (command, SZ_LINE, TY_CHAR)
+
+	while (clgcur ("cursor", xx, yy, wcs, key, Memc[command],
+	    SZ_LINE) != EOF) {
+
+	    x = double(xx)
+	    y = double(yy)
+
+	    # Switch on command, take appropriate action.
+	    switch (key) {
+	    case '?':
+	        # List options.
+		call gpagefile (PDM_GP(pdmp), HELP, PROMPT)
+	    case ':':
+		# Colon command.
+		call pdm_colon (pdmp, Memc[command], ptype, infile,
+		    period, flip)
+	    case 'h':
+		# Graph the data.
+		call pdm_dplot (pdmp, infile, flip)
+		ptype = DATAPLOT
+	    case 'f':
+		# Call icfit on the data.
+		if (ptype == DATAPLOT) {
+		    # Set the min/max ordinate values for icfit.
+		    call alimd (PDM_X(pdmp,1), PDM_NPT(pdmp), xmin, xmax)
+		    call ic_putr (PDM_ICD(pdmp), "xmin", real(xmin))
+		    call ic_putr (PDM_ICD(pdmp), "xmax", real(xmax))
+
+		    # Allocate a temporary weights array and fill it with the
+		    # in-use array values.
+
+		    call smark (sptemp)
+		    call salloc (weights, PDM_NPT(pdmp), TY_DOUBLE)
+		    do i = 1, PDM_NPT(pdmp)
+			Memd[weights+i-1] = double(PDM_INUSE(pdmp,i))
+
+		    # Call icfit.
+		    if (PDM_NPT(pdmp) >= 2)
+	    		call icg_fitd (PDM_ICD(pdmp), PDM_GP(pdmp), "cursor",
+			    PDM_GT(pdmp), PDM_CVD(pdmp), PDM_X(pdmp,1),
+			    PDM_DY(pdmp,1), Memd[weights], PDM_NPT(pdmp))
+
+		    # Recover the weights array back into the in-use array.
+		    do i = 1, PDM_NPT(pdmp)
+			 PDM_INUSE(pdmp,i) = int(Memd[weights+i-1])
+
+		    call sfree (sptemp)
+
+		    # Replot.
+		    call pdm_dplot (pdmp, infile, flip)
+		    call sfree (sptemp)
+		} else if (ptype == PHASEPLOT) {
+		    # Call icfit on the phases.
+		    call pdm_fitphase (pdmp)
+
+		    # Replot.
+		    call pdm_pplot (pdmp, period, infile, flip)
+		} else
+		    call printf ("Can't fit a THETA plot \007\n")
+	    case 'i':
+		# Theta vs. frequency plot.
+		if (PDM_PMIN(pdmp) < EPSILOND && PDM_PMAX(pdmp) < EPSILOND)
+		    call pdm_minmaxp(pdmp)
+
+		# Calculate theta.
+		call pdm_alltheta (pdmp, THETAFPLOT)
+
+		# Graph theta vs frequency.
+		call pdm_tplot (pdmp, THETAFPLOT, infile)
+		ptype = THETAFPLOT
+	    case 'k':
+		# Theta vs. period plot.
+		if (PDM_PMIN(pdmp) < EPSILOND && PDM_PMAX(pdmp) < EPSILOND)
+		    call pdm_minmaxp(pdmp)
+
+		# Calculate theta.
+		call pdm_alltheta (pdmp, THETAPPLOT)
+
+		# Graph theta vs frequency.
+		call pdm_tplot (pdmp, THETAPPLOT, infile)
+		ptype = THETAPPLOT
+	    case 'p':
+		# Phase curve plot.
+		if (ptype == THETAPPLOT) {
+		    # Find the epoch, calculate the phases w.r.t. this epoch.
+		    call pdm_ampep (pdmp, x)
+		    call pdm_phase(pdmp, x, PDM_EPOCH(pdmp))
+		    call pdm_pplot (pdmp, x, infile, flip)
+		    ptype = PHASEPLOT
+		    period = x
+		} else if (ptype == THETAFPLOT) {
+		    # Find the epoch, calculate the phases w.r.t. this epoch.
+		    call pdm_ampep (pdmp, 1.0d+0/x)
+		    call pdm_phase(pdmp, 1.0d+0/x, PDM_EPOCH(pdmp))
+		    call pdm_pplot (pdmp, 1.0d+0/x, infile, flip)
+		    ptype = PHASEPLOT
+		    period = 1.0d+0/x
+		} else
+		    call printf ("Wrong type of plot on screen for p key\007\n")
+	    case 'd':
+		# Delete the point nearest the cursor.
+		index = pdm_delete (pdmp, x, y, ptype)
+	    case 'u':
+		# Undelete the point nearest the cursor.
+		index = pdm_undelete (pdmp, x, y, ptype)
+	    case 'j':
+		# Subtract the fit from the data and use residuals.
+		if (PDM_CVD(pdmp) == NULL)
+		    call printf ("Fit has not been done. \007\n")
+		else if (PDM_RESID(pdmp) == YES)
+		    call printf ("Already using residuals. \007\n")
+		else {
+		    # For each point, calculate the fit function and subtract
+		    # it from the data.
+
+		    do i = 1, PDM_NPT(pdmp) {
+			PDM_DY(pdmp,i) = PDM_DY(pdmp,i) -
+			    dcveval (PDM_CVD(pdmp), PDM_X(pdmp,i))
+		    }
+		    PDM_RESID(pdmp) = YES
+		    if (ptype == DATAPLOT)
+		        call pdm_dplot (pdmp, infile, flip)
+		}
+	    case 's':
+		# Set sample regions with the cursor.
+		if (ptype == DATAPLOT) {
+		    if (stridxs ("*", PDM_SAMPLE(pdmp)) > 0)
+		        PDM_SAMPLE(pdmp) = EOS
+
+		    rx1 = x
+		    call printf ("again:\n")
+		    if (clgcur ("cursor", xx, yy, wcs, key, Memc[command],
+			SZ_LINE) == EOF)
+		        break
+		    rx2 = double(xx)
+
+		    call sprintf (Memc[command], SZ_LINE, " %g:%g")
+		        call pargd (rx1)
+		        call pargd (rx2)
+
+		    call strcat (Memc[command], PDM_SAMPLE(pdmp), SZ_LINE)
+		    call rg_gxmarkd (PDM_GP(pdmp), PDM_SAMPLE(pdmp),
+			PDM_X(pdmp,1), PDM_NPT(pdmp), 1)
+		    call printf (" \n")
+		    call gflush (PDM_GP(pdmp))
+		} else
+		    call printf ("Wrong type of plot for s key \007\n")
+	    case 't':
+		# Initialize the sample string and erase from the graph.
+		if (ptype == DATAPLOT) {
+		    call rg_gxmarkd (PDM_GP(pdmp), PDM_SAMPLE(pdmp),
+			PDM_X(pdmp,1), PDM_NPT(pdmp), 0)
+		}
+		call gflush (PDM_GP(pdmp))
+		call strcpy ("*", PDM_SAMPLE(pdmp), SZ_LINE)
+		call gflush (PDM_GP(pdmp))
+	    case 'g':
+		# Significance of theta at cursor x position.
+		if (ptype == THETAPPLOT) {
+	    	    # Calculate significance at cursor x position (per).
+		        signif = pdm_signif (pdmp, x)
+	    	    # Print at bottom of screen.
+		    call printf ("Significance at cursor = %g\n")
+			call pargd (signif)
+		} else if (ptype == THETAFPLOT) {
+	    	    # Calculate significance at cursor x position (per).
+		        signif = pdm_signif (pdmp, 1.0d+0/x)
+		    # Print at bottom of screen.
+		    call printf ("Significance at cursor = %g\n")
+			call pargd (signif)
+		} else if (ptype == PHASEPLOT) {
+		    # Calculate significance at current period/frequency
+		        signif = pdm_signif (pdmp, period)
+		    # Print at bottom of screen
+		    call printf ("Significance of this period = %g\n")
+			call pargd (signif)
+		} else {  
+		    # Data plot.
+		    call printf ("Wrong type of plot for g key \007\n")
+		}
+	    case 'a':
+		# Amplitude and epoch at cursor x position.
+		if (ptype == THETAPPLOT) {
+		    # Calculate ampl & epoch at cursor x position.
+		    call pdm_ampep (pdmp, x)
+		    # Print at bottom of screen.
+		    call printf ("amplitude of data = %g, epoch = %g\n")
+			call pargd (PDM_AMPL(pdmp))
+			call pargd (PDM_EPOCH(pdmp))
+		} else if (ptype == THETAFPLOT) {
+		    # Calculate ampl & epoch at cursor x position.
+		    call pdm_ampep (pdmp, 1.0d+0/x)
+		    # Print at bottom of screen.
+		    call printf ("amplitude of data = %g, epoch = %g\n")
+			call pargd (PDM_AMPL(pdmp))
+			call pargd (PDM_EPOCH(pdmp))
+		} else if (ptype == PHASEPLOT) {
+		    # Calculate ampl & epoch at current period/frequency
+		    call pdm_ampep (pdmp, period)
+		    # Print at bottom of screen.
+		    call printf ("amplitude of data = %g, epoch = %g\n")
+			call pargd (PDM_AMPL(pdmp))
+			call pargd (PDM_EPOCH(pdmp))
+		} else {
+		    # Data plot.
+		    call printf ("Wrong type of plot for a key \007\n")
+		}
+	    case ',':
+		# Set minp or minf to cursor x position.
+		if (ptype == THETAPPLOT) {
+		    PDM_PMIN(pdmp) = x
+		    PDM_FMIN(pdmp) = 1.0d+0/x
+		    # Print at bottom of screen.
+		    call printf ("minp now %g\n")
+			call pargd (PDM_PMIN(pdmp))
+		} else if (ptype == THETAFPLOT) {
+		    PDM_FMAX(pdmp) = x
+		    PDM_PMAX(pdmp) = 1.0d+0/x
+		    # Print at bottom of screen.
+		    call printf ("minf now %g\n")
+			call pargd (PDM_FMAX(pdmp))
+		} else {
+		    # Data plot or phase plot.
+		    call printf ("Wrong type of plot for , key \007\n")
+		}
+	    case '.':
+		# Set maxp or maxf to cursor x position.
+		if (ptype == THETAPPLOT) {
+		    PDM_PMAX(pdmp) = x
+		    PDM_FMAX(pdmp) = 1.0d+0/x
+		    # Print at bottom of screen.
+		    call printf ("maxp now %g\n")
+			call pargd (PDM_PMAX(pdmp))
+		} else if (ptype == THETAFPLOT) {
+		    PDM_FMIN(pdmp) = x
+		    PDM_PMIN(pdmp) = 1.0d+0/x
+		    # Print at bottom of screen.
+		    call printf ("maxf now %g\n")
+			call pargd (PDM_FMIN(pdmp))
+		} else {
+		    # Data plot or phase plot.
+		    call printf ("Wrong type of plot for . key \007\n")
+		}
+	    case 'm':
+		# Mark range and find minimum in this range.
+		if (ptype == THETAFPLOT || ptype == THETAPPLOT) {
+		    rx1 = x
+		    call printf ("again:\n")
+		    if (clgcur ("cursor", xx, yy, wcs, key, Memc[command],
+			SZ_LINE) == EOF)
+		        break
+		    rx2 = double(xx)
+		    index = pdm_findmin(pdmp, ptype, rx1, rx2, 1,
+			PDM_NTHPT(pdmp))
+	    	    PDM_MINR(pdmp) = PDM_XTH(pdmp,index)
+		    call printf ("period at minimum = %g, frequency = %g\n")
+			call pargd (PDM_XTH(pdmp,index))
+			call pargd (1.0d+0/PDM_XTH(pdmp,index))
+
+		} else
+		    call printf ("Wrong type of plot for m key. \007\n")
+	    case 'r':
+		# Check type of plot and replot.
+		if (ptype == DATAPLOT)
+		    call pdm_dplot (pdmp, infile, flip)
+		if (ptype == THETAPPLOT)
+		    call pdm_tplot (pdmp, THETAPPLOT, infile)
+		if (ptype == THETAFPLOT)
+		    call pdm_tplot (pdmp, THETAFPLOT, infile)
+		if (ptype == PHASEPLOT)
+		    call pdm_pplot (pdmp, period, infile, flip)
+	    case 'e':
+		# Toggle error bars.
+		if (PDM_EB(pdmp) == NO)
+		    PDM_EB(pdmp) = YES
+		else
+		    PDM_EB(pdmp) = NO
+
+		# Check type of plot and replot.
+		if (ptype == PHASEPLOT)
+		    call pdm_pplot (pdmp, period, infile, flip)
+	    case 'x':
+		# Remove a trend from the data by fitting a straight line to
+		# the data and removing it.
+
+		# Set the min/max ordinate values for icfit.
+		call alimd (PDM_X(pdmp,1), PDM_NPT(pdmp), xmin, xmax)
+		call dcvinit (PDM_CVD(pdmp), SPLINE1, 1, xmin, xmax)
+
+		# Allocate a weights array and fill it.
+		call smark (sptemp)
+		call salloc (weights, PDM_NPT(pdmp), TY_DOUBLE)
+		do i = 1, PDM_NPT(pdmp)
+		    Memd[weights+i-1] = PDM_INUSE(pdmp,i)
+
+		call dcvfit (PDM_CVD(pdmp), PDM_X(pdmp,1), PDM_DY(pdmp,1),
+		    Memd[weights], PDM_NPT(pdmp), WTS_USER, ier)
+		if (ier != 0) {
+		    call eprintf ("error in dcvfit\n")
+		    call erract (EA_WARN)
+		}
+
+		# Subtract the fit from the data.
+		if (PDM_RESID(pdmp) == YES) {
+		    call printf ("Already using residuals. \007\n")
+		    next
+		} else {
+		    # For each point, calculate the fit function and subtract
+		    # it from the data.
+
+		    do i = 1, PDM_NPT(pdmp) {
+			PDM_DY(pdmp,i) = PDM_DY(pdmp,i) -
+			    dcveval (PDM_CVD(pdmp), PDM_X(pdmp,i))
+		    }
+		    PDM_RESID(pdmp) = YES
+		}
+
+		call dcvfree (PDM_CVD(pdmp))
+		call sfree (sptemp)
+
+		# Replot.
+		call pdm_dplot (pdmp, infile, flip)
+		ptype = DATAPLOT
+	    case 'z':
+		# Flip the y-axis scale.
+		flip = !flip
+	    case 'q':
+		# Quit. Exit PDM.
+		if (gqverify() == YES)
+		    break
+	    default:
+		# Error: unknown command: ring bell.
+		call printf ("\007\n")
+		call printf ("? for help or (h,f,i,k,p")
+		call printf (",d,u,j,s,t,g,a,m,r,x,q,:)\n")
+	    }
+	}
+
+	call flush (STDOUT)
+	call sfree (sp)
+end