aboutsummaryrefslogtreecommitdiff
path: root/noao/digiphot/daophot/daolib/dpgetapert.x
diff options
context:
space:
mode:
Diffstat (limited to 'noao/digiphot/daophot/daolib/dpgetapert.x')
-rw-r--r--noao/digiphot/daophot/daolib/dpgetapert.x530
1 files changed, 530 insertions, 0 deletions
diff --git a/noao/digiphot/daophot/daolib/dpgetapert.x b/noao/digiphot/daophot/daolib/dpgetapert.x
new file mode 100644
index 00000000..8268373f
--- /dev/null
+++ b/noao/digiphot/daophot/daolib/dpgetapert.x
@@ -0,0 +1,530 @@
+include <tbset.h>
+include "../../lib/ptkeysdef.h"
+include "../lib/daophotdef.h"
+include "../lib/apseldef.h"
+
+# DP_WGETAPERT -- Read the aperture photometry results and transform the
+# input coordinates to the appropriate coordinate system. Works with
+# either the "old" APPHOT files or the new ST Tables.
+
+procedure dp_wgetapert (dao, im, apd, max_nstars, old_ap)
+
+pointer dao # pointer to the DAOPHOT structure
+pointer im # the input image descriptor
+int apd # input photometry file descriptor
+int max_nstars # maximum number of stars
+bool old_ap # YES indicates old APPHOT file
+
+pointer apsel
+int dp_stati()
+
+begin
+ # Get the stars.
+ call dp_getapert (dao, apd, max_nstars, old_ap)
+
+ # Transform the coordinates if necessary.
+ apsel = DP_APSEL (dao)
+ if (dp_stati (dao, WCSIN) != WCS_LOGICAL)
+ call dp_win (dao, im, Memr[DP_APXCEN(apsel)],
+ Memr[DP_APYCEN(apsel)], Memr[DP_APXCEN(apsel)],
+ Memr[DP_APYCEN(apsel)], DP_APNUM(apsel))
+end
+
+
+# DP_GETAPERT -- Read the aperture photometry results. Works with
+# either the "old" APPHOT files or the new ST Tables.
+
+procedure dp_getapert (dao, apd, max_nstars, old_ap)
+
+pointer dao # pointer to the DAOPHOT structure
+int apd # input Photometry file descriptor
+int max_nstars # maximum number of stars
+bool old_ap # YES indicates old APPHOT file
+
+int nstars
+pointer apsel
+int tbpsta(), dp_goldap(), dp_gtabphot()
+
+begin
+ # Get APSEL pointer.
+ apsel = DP_APSEL (dao)
+
+ # Get the required memory.
+ Memi[DP_APRESULT(apsel)] = DP_PAPID
+ Memi[DP_APRESULT(apsel)+1] = DP_PAPXCEN
+ Memi[DP_APRESULT(apsel)+2] = DP_PAPYCEN
+ Memi[DP_APRESULT(apsel)+3] = DP_PAPMAG1
+ Memi[DP_APRESULT(apsel)+4] = DP_PAPSKY
+
+ # Get the results.
+ if (old_ap) {
+ call dp_memapsel (dao, Memi[DP_APRESULT(apsel)], NAPRESULT,
+ max_nstars)
+ nstars = dp_goldap (apd, dao, max_nstars)
+ } else {
+ call dp_memapsel (dao, Memi[DP_APRESULT(apsel)], NAPRESULT,
+ tbpsta (apd, TBL_NROWS))
+ nstars = dp_gtabphot (apd, dao, tbpsta (apd, TBL_NROWS))
+ }
+
+ # Reallocate to save space if appropopriate.
+ if (nstars < max_nstars)
+ call dp_rmemapsel (dao, Memi[DP_APRESULT(apsel)], NAPRESULT,
+ nstars)
+end
+
+
+# DP_MEMAPSEL -- Procedure to allocate memory for the apselect strucuture.
+
+procedure dp_memapsel (dao, fields, nfields, max_nstars)
+
+pointer dao # pointer to the daophot strucuture
+int fields[ARB] # array of fields
+int nfields # number of fields to allocate space for
+int max_nstars # maximum number of stars
+
+int i
+pointer apsel
+
+begin
+ apsel = DP_APSEL(dao)
+
+ # Allocate space for results.
+ do i = 1, nfields {
+ switch (fields[i]) {
+ case DP_PAPID:
+ if (DP_APID(apsel) != NULL)
+ call mfree (DP_APID(apsel), TY_INT)
+ call malloc (DP_APID(apsel), max_nstars, TY_INT)
+
+ case DP_PAPXCEN:
+ if (DP_APXCEN(apsel) != NULL)
+ call mfree (DP_APXCEN(apsel), TY_REAL)
+ call malloc (DP_APXCEN(apsel), max_nstars, TY_REAL)
+
+ case DP_PAPYCEN:
+ if (DP_APYCEN(apsel) != NULL)
+ call mfree (DP_APYCEN(apsel), TY_REAL)
+ call malloc (DP_APYCEN(apsel), max_nstars, TY_REAL)
+
+ case DP_PAPSKY:
+ if (DP_APMSKY(apsel) != NULL)
+ call mfree (DP_APMSKY(apsel), TY_REAL)
+ call malloc (DP_APMSKY(apsel), max_nstars, TY_REAL)
+
+ case DP_PAPMAG1:
+ if (DP_APMAG(apsel) != NULL)
+ call mfree (DP_APMAG(apsel), TY_REAL)
+ call malloc (DP_APMAG(apsel), max_nstars, TY_REAL)
+
+ case DP_PAPGROUP:
+ if (DP_APGROUP(apsel) != NULL)
+ call mfree (DP_APGROUP(apsel), TY_INT)
+ #call malloc (DP_APGROUP(apsel), max_nstars, TY_INT)
+
+ case DP_PAPMERR1:
+ if (DP_APERR(apsel) != NULL)
+ call mfree (DP_APERR(apsel), TY_REAL)
+ call malloc (DP_APERR(apsel), max_nstars, TY_REAL)
+
+ case DP_PAPNITER:
+ if (DP_APNITER(apsel) != NULL)
+ call mfree (DP_APNITER(apsel), TY_INT)
+ #call malloc (DP_APNITER(apsel), max_nstars, TY_INT)
+
+ case DP_PAPCHI:
+ if (DP_APCHI(apsel) != NULL)
+ call mfree (DP_APCHI(apsel), TY_REAL)
+ call malloc (DP_APCHI(apsel), max_nstars, TY_REAL)
+
+ case DP_PAPSHARP:
+ if (DP_APSHARP(apsel) != NULL)
+ call mfree (DP_APSHARP(apsel), TY_REAL)
+ call malloc (DP_APSHARP(apsel), max_nstars, TY_REAL)
+ }
+ }
+end
+
+
+# DP_GOLDAP -- Read in the photometry from an old style APPHOT file
+
+int procedure dp_goldap (apd, dao, max_nstars)
+
+int apd # the input file descriptor
+pointer dao # pointer to the daophot structure
+int max_nstars # maximum number of stars
+
+int nstars, bufsize, stat
+pointer apsel, apkey, sp, fields
+int dp_apsel()
+
+begin
+ # Define the point to the apselect structure.
+ apsel = DP_APSEL (dao)
+
+ # Allocate some temporary space.
+ call smark (sp)
+ call salloc (fields, SZ_LINE, TY_CHAR)
+
+ # Initialize the keyword structure.
+ call pt_kyinit (apkey)
+
+ # Set up the fields to be retrieved.
+ call dp_gappsf (Memi[DP_APRESULT(apsel)], Memc[fields], NAPRESULT)
+
+ # Now read in the results.
+ nstars = 0
+ bufsize = max_nstars
+ repeat {
+
+ # Read in a group of stars.
+ while (nstars < bufsize) {
+ stat = dp_apsel (apkey, apd, Memc[fields],
+ Memi[DP_APRESULT(apsel)], Memi[DP_APID(apsel)+nstars],
+ Memr[DP_APXCEN(apsel)+nstars],
+ Memr[DP_APYCEN(apsel)+nstars],
+ Memr[DP_APMSKY(apsel)+nstars],
+ Memr[DP_APMAG(apsel)+nstars])
+ if (stat == EOF)
+ break
+ nstars = nstars + 1
+ }
+
+ # Check the buffer size.
+ if (stat == EOF)
+ break
+ bufsize = bufsize + max_nstars
+ call dp_rmemapsel (dao, Memi[DP_APRESULT(apsel)], NAPRESULT,
+ bufsize)
+
+ }
+ DP_APNUM(apsel) = nstars
+
+ # Free the keyword structure.
+ call pt_kyfree (apkey)
+ call sfree (sp)
+
+ return (nstars)
+end
+
+
+# DP_GTABPHOT -- Read in the complete photometry from an ST table.
+
+int procedure dp_gtabphot (tp, dao, max_nstars)
+
+pointer tp # table descriptor
+pointer dao # pointer to daophot structure
+int max_nstars # maximum number of stars
+
+bool nullflag
+int record, index, nrow
+pointer apsel, idpt, xcenpt, ycenpt, magpt, skypt
+int tbpsta()
+
+begin
+ # Define the point to the apselect structure.
+ apsel = DP_APSEL (dao)
+
+ # Find the column pointers
+ call tbcfnd (tp, ID, idpt, 1)
+ if (idpt == NULL)
+ call tbcfnd (tp, "ID", idpt, 1)
+ if (idpt == NULL)
+ call printf ("Error reading ID.\n")
+
+ call tbcfnd (tp, XCENTER, xcenpt, 1)
+ if (xcenpt == NULL)
+ call tbcfnd (tp, "XCENTER", xcenpt, 1)
+ if (xcenpt == NULL)
+ call printf ("Error reading XCENTER.\n")
+
+ call tbcfnd (tp, YCENTER, ycenpt, 1)
+ if (ycenpt == NULL)
+ call tbcfnd (tp, "YCENTER", ycenpt, 1)
+ if (ycenpt == NULL)
+ call printf ("Error reading YCENTER.\n")
+
+ call tbcfnd (tp, MAG, magpt, 1)
+ if (magpt == NULL)
+ call tbcfnd (tp, APMAG, magpt, 1)
+ if (magpt == NULL)
+ call printf ("Error reading MAG.\n")
+
+ call tbcfnd (tp, SKY, skypt, 1)
+ if (skypt == NULL)
+ call tbcfnd (tp, SKY, skypt, 1)
+ if (skypt == NULL)
+ call printf ("Error reading SKY.\n")
+
+
+ # Get the results ignoring any record with ID = NULL.
+ nrow = min (tbpsta (tp, TBL_NROWS), max_nstars)
+ index = 0
+ do record = 1, nrow {
+
+ # Check the ID record.
+ call tbrgti (tp, idpt, Memi[DP_APID(apsel)+index], nullflag, 1,
+ record)
+ if (nullflag)
+ next
+
+ # Read the remaining records.
+ call tbrgtr (tp, xcenpt, Memr[DP_APXCEN(apsel)+index], nullflag, 1,
+ record)
+ call tbrgtr (tp, ycenpt, Memr[DP_APYCEN(apsel)+index], nullflag, 1,
+ record)
+ call tbrgtr (tp, magpt, Memr[DP_APMAG(apsel)+index], nullflag, 1,
+ record)
+ call tbrgtr (tp, skypt, Memr[DP_APMSKY(apsel)+index], nullflag, 1,
+ record)
+ index = index + 1
+ }
+
+ DP_APNUM(apsel) = index
+
+ return (index)
+end
+
+
+# DP_RMEMAPSEL -- Procedure to reallocate memory for the apselect strucuture.
+
+procedure dp_rmemapsel (dao, fields, nfields, max_nstars)
+
+pointer dao # pointer to the daophot strucuture
+int fields[ARB] # integer fields
+int nfields # number of fields
+int max_nstars # maximum number of stars
+
+int i
+pointer apsel
+
+begin
+ # Reallocate space for results.
+ apsel = DP_APSEL(dao)
+
+ do i = 1, nfields {
+ switch (fields[i]) {
+ case DP_PAPID:
+ call realloc (DP_APID(apsel), max_nstars, TY_INT)
+ case DP_PAPXCEN:
+ call realloc (DP_APXCEN(apsel), max_nstars, TY_REAL)
+ case DP_PAPYCEN:
+ call realloc (DP_APYCEN(apsel), max_nstars, TY_REAL)
+ case DP_PAPSKY:
+ call realloc (DP_APMSKY(apsel), max_nstars, TY_REAL)
+ case DP_PAPGROUP:
+ #call realloc (DP_APGROUP(apsel), max_nstars, TY_INT)
+ case DP_PAPMAG1:
+ call realloc (DP_APMAG(apsel), max_nstars, TY_REAL)
+ case DP_PAPMERR1:
+ call realloc (DP_APERR(apsel), max_nstars, TY_REAL)
+ case DP_PAPNITER:
+ #call realloc (DP_APNITER(apsel), max_nstars, TY_INT)
+ case DP_PAPSHARP:
+ call realloc (DP_APSHARP(apsel), max_nstars, TY_REAL)
+ case DP_PAPCHI:
+ call realloc (DP_APCHI(apsel), max_nstars, TY_REAL)
+ }
+ }
+end
+
+
+# DP_GAPPSF -- Set up the structures necessary for retrieving the
+# aperture phtometery results needed for the PSF task.
+
+procedure dp_gappsf (fields, sel_fields, max_nfields)
+
+int fields[ARB] # array of selected fields
+char sel_fields[ARB] # names of selected containing fields
+int max_nfields # maximum number of fields selected
+
+int i
+int strlen()
+
+begin
+ sel_fields[1] = EOS
+
+ do i = 1, max_nfields {
+ switch (fields[i]) {
+ case DP_PAPID:
+ call sprintf (sel_fields[strlen(sel_fields)+1], SZ_LINE, "%s ")
+ call pargstr (ID)
+ case DP_PAPXCEN:
+ call sprintf (sel_fields[strlen(sel_fields)+1], SZ_LINE, "%s ")
+ call pargstr (XCENTER)
+ case DP_PAPYCEN:
+ call sprintf (sel_fields[strlen(sel_fields)+1], SZ_LINE, "%s ")
+ call pargstr (YCENTER)
+ case DP_PAPSKY:
+ call sprintf (sel_fields[strlen(sel_fields)+1], SZ_LINE, "%s ")
+ call pargstr (SKY)
+ case DP_PAPMAG1:
+ call sprintf (sel_fields[strlen(sel_fields)+1], SZ_LINE, "%s ")
+ call pargstr (APMAG)
+ }
+ }
+
+ if (sel_fields[1] != EOS)
+ sel_fields[strlen(sel_fields)] = EOS
+end
+
+
+# DP_APSEL -- Select records from an apphot/daophot text file.
+
+int procedure dp_apsel (key, fd, fields, indices, id, x, y, sky, mag)
+
+pointer key # pointer to key structure
+int fd # text file descriptor
+char fields[ARB] # fields to be output
+int indices[ARB] # indices of fields
+int id # star id number
+real x # x center
+real y # y center
+real sky # sky value
+real mag # magnitude
+
+int nchars, nunique, uunique, funique, ncontinue, recptr
+int first_rec, nselect, record
+pointer line
+int getline(), strncmp(), pt_choose()
+
+data first_rec /YES/
+
+begin
+ # Initialize the file read.
+ if (first_rec == YES) {
+ nunique = 0
+ uunique = 0
+ funique = 0
+ nselect = 0
+ record = 0
+ call malloc (line, SZ_LINE, TY_CHAR)
+ }
+
+ ncontinue = 0
+ recptr = 1
+
+ # Loop over the text file records.
+ repeat {
+
+ # Read in a line of the text file.
+ nchars = getline (fd, Memc[line])
+ if (nchars == EOF)
+ break
+
+ # Determine the type of record.
+ if (Memc[line] == KY_CHAR_POUND) {
+
+ if (strncmp (Memc[line], KY_CHAR_KEYWORD, KY_LEN_STR) == 0) {
+ call pt_kyadd (key, Memc[line], nchars)
+ } else if (strncmp (Memc[line], KY_CHAR_NAME,
+ KY_LEN_STR) == 0) {
+ nunique = nunique + 1
+ call pt_kname (key, Memc[line], nchars, nunique)
+ } else if (strncmp (Memc[line], KY_CHAR_UNITS,
+ KY_LEN_STR) == 0) {
+ uunique = uunique + 1
+ call pt_knunits (key, Memc[line], nchars, uunique)
+ } else if (strncmp (Memc[line], KY_CHAR_FORMAT,
+ KY_LEN_STR) == 0) {
+ funique = funique + 1
+ call pt_knformats (key, Memc[line], nchars, funique)
+ }
+ } else if (Memc[line] == KY_CHAR_NEWLINE) {
+ # skip blank lines
+
+ } else {
+
+ # Construct the table record.
+ call pt_mkrec (key, Memc[line], nchars, first_rec, recptr,
+ ncontinue)
+
+ # Construct output record when there is no continuation char.
+ if (Memc[line+nchars-2] != KY_CHAR_CONT) {
+
+ # Select the appropriate records.
+ if (nselect <= 0) {
+ nselect = pt_choose (key, fields)
+ if (nselect <= 0)
+ break
+ }
+
+ # Construct the output record by moving selected fields
+ # into data structures.
+ call dp_getap (key, indices, id, x, y, sky, mag)
+ record = record + 1
+ first_rec = NO
+
+ # Record is complete so exist the loop.
+ break
+ }
+ }
+
+ }
+
+ if (nchars == EOF || nselect <= 0) {
+ first_rec = YES
+ nunique = 0
+ uunique = 0
+ funique = 0
+ nselect = 0
+ call mfree (line, TY_CHAR)
+ return (EOF)
+ } else
+ return (record)
+end
+
+
+# DP_GETAP -- Decode the selected daophot/apphot values.
+
+procedure dp_getap (key, indices, id, x, y, sky, mag)
+
+pointer key # pointer to keys strucuture
+int indices[ARB] # index array
+int id # star id
+real x # x position
+real y # y position
+real sky # sky value
+real mag # magnitude
+
+int i, index, elem, maxch, kip, ip
+int ctoi(), ctor()
+char buffer[SZ_LINE]
+
+begin
+ do i = 1, KY_NSELECT(key) {
+
+ # Find the key.
+ index = Memi[KY_SELECT(key)+i-1]
+ elem = Memi[KY_ELEM_SELECT(key)+i-1]
+ maxch = Memi[KY_LEN_SELECT(key)+i-1]
+ kip = Memi[KY_PTRS(key)+index-1] + (elem - 1) * maxch
+
+ # Extract the appropriate field.
+ call amovc (Memc[kip], buffer, maxch)
+ buffer[maxch+1] = EOS
+
+ # Decode the output value.
+ ip = 1
+ switch (indices[i]) {
+ case DP_PAPID:
+ if (ctoi (buffer, ip, id) <= 0)
+ id = 0
+ case DP_PAPXCEN:
+ if (ctor (buffer, ip, x) <= 0)
+ x = INDEFR
+ case DP_PAPYCEN:
+ if (ctor (buffer, ip, y) <= 0)
+ y = INDEFR
+ case DP_PAPSKY:
+ if (ctor (buffer, ip, sky) <= 0)
+ sky = INDEFR
+ case DP_PAPMAG1:
+ if (ctor (buffer, ip, mag) <= 0)
+ mag = INDEFR
+ default:
+ call printf ("Error reading the APPHOT results.\n")
+ }
+
+ }
+end
generated by cgit (git 2.34.1) at 2025-09-20 08:59:51 -0400