aboutsummaryrefslogtreecommitdiff
path: root/noao/astutil/pdm/pdmcursor.x
diff options
context:
space:
mode:
Diffstat (limited to 'noao/astutil/pdm/pdmcursor.x')
-rw-r--r--noao/astutil/pdm/pdmcursor.x383
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
generated by cgit (git 2.34.1) at 2025-09-20 05:23:34 -0400