From fa080de7afc95aa1c19a6e6fc0e0708ced2eadc4 Mon Sep 17 00:00:00 2001 From: Joseph Hunkeler Date: Wed, 8 Jul 2015 20:46:52 -0400 Subject: Initial commit --- noao/digiphot/apphot/polyphot/apyget.x | 306 +++++++++++++++++++++++++++++++++ 1 file changed, 306 insertions(+) create mode 100644 noao/digiphot/apphot/polyphot/apyget.x (limited to 'noao/digiphot/apphot/polyphot/apyget.x') diff --git a/noao/digiphot/apphot/polyphot/apyget.x b/noao/digiphot/apphot/polyphot/apyget.x new file mode 100644 index 00000000..a990b6f8 --- /dev/null +++ b/noao/digiphot/apphot/polyphot/apyget.x @@ -0,0 +1,306 @@ +include +include "../lib/apphot.h" +include "../lib/fitsky.h" +include "../lib/polyphot.h" + +# AP_YGET -- Procedure to fetch the coordinates of the vertices of a +# polygon from a text file. + +int procedure ap_yget (py, im, fd, delim, x, y, max_nvertices) + +pointer py # polyphot structure +pointer im # the input image descriptor +int fd # polygon file descriptor +int delim # delimiter character +real x[ARB] # x coords of vertices +real y[ARB] # y coords of vertices +int max_nvertices # maximum number of vertices + +real xc, yc, r2max, r2, xtemp, ytemp +pointer sp, str +int i, nvertices, nreal, stat +char marker +real asumr() +int fscan(), nscan(), strncmp(), apstati() + +begin + call smark (sp) + call salloc (str, SZ_LINE, TY_CHAR) + + # Print prompts if input file is STDIN. + call fstats (fd, F_FILENAME, Memc[str], SZ_LINE) + if (strncmp ("STDIN", Memc[str], 5) == 0) { + call printf ("Type x and y coordinates of vertex\n") + call printf ("\t1 vertex per line\n") + call printf ("\t; to end polygon\n") + call printf ("\t^Z to end list\n") + call flush (STDOUT) + } + + # Get the polygon. + nvertices = 0 + stat = fscan (fd) + while (stat != EOF) { + + # Read the vertices from the file + call gargr (xtemp) + call gargr (ytemp) + nreal = nscan () + if (nreal != 2) { + call reset_scan () + call gargc (marker) + if (int (marker) == delim) + break + nreal = nscan () + } + + # Store the vertices. + if (nreal >= 2) { + nvertices = nvertices + 1 + if (nvertices <= max_nvertices) { + x[nvertices] = xtemp + y[nvertices] = ytemp + } + } + + stat = fscan (fd) + } + + if (nvertices == 0 && stat == EOF) { + + nvertices = EOF + call apsetr (py, PYXMEAN, INDEFR) + call apsetr (py, PYYMEAN, INDEFR) + call apsetr (py, OPYXMEAN, INDEFR) + call apsetr (py, OPYYMEAN, INDEFR) + call apsetr (py, PYCX, INDEFR) + call apsetr (py, PYCY, INDEFR) + call apsetr (py, OPYCX, INDEFR) + call apsetr (py, OPYCY, INDEFR) + call apsetr (py, PYMINRAD, INDEFR) + call apseti (py, PYNVER, 0) + + } else if (nvertices <= 2) { + + call apsetr (py, PYXMEAN, INDEFR) + call apsetr (py, PYYMEAN, INDEFR) + call apsetr (py, OPYXMEAN, INDEFR) + call apsetr (py, OPYYMEAN, INDEFR) + call apsetr (py, PYCX, INDEFR) + call apsetr (py, PYCY, INDEFR) + call apsetr (py, OPYCX, INDEFR) + call apsetr (py, OPYCY, INDEFR) + call apsetr (py, PYMINRAD, INDEFR) + call apseti (py, PYNVER, 0) + nvertices = 0 + + } else { + + # Add in the last vertex. + x[nvertices+1] = x[1] + y[nvertices+1] = y[1] + + # Transform the input coordinates. + switch (apstati(py,WCSIN)) { + case WCS_WORLD, WCS_PHYSICAL: + call ap_itol (py, x, y, x, y, nvertices + 1) + case WCS_TV: + call ap_vtol (im, x, y, x, y, nvertices + 1) + default: + ; + } + + # Compute the mean polygon coordinates. + xc = asumr (x, nvertices) / nvertices + yc = asumr (y, nvertices) / nvertices + call apsetr (py, PYXMEAN, xc) + call apsetr (py, PYYMEAN, yc) + call apsetr (py, PYCX, xc) + call apsetr (py, PYCY, yc) + + # Set the minimum size of the sky annulus. + r2max = 0.0 + do i = 1, nvertices { + r2 = (x[i] - xc) ** 2 + (y[i] - yc) ** 2 + if (r2 > r2max) + r2max = r2 + } + + switch (apstati(py,WCSOUT)) { + case WCS_WORLD, WCS_PHYSICAL: + call ap_ltoo (py, xc, yc, xc, yc, 1) + case WCS_TV: + call ap_ltov (im, xc, yc, xc, yc, 1) + default: + ; + } + call apsetr (py, OPYXMEAN, xc) + call apsetr (py, OPYYMEAN, yc) + call apsetr (py, OPYCX, xc) + call apsetr (py, OPYCY, yc) + + call apsetr (py, PYMINRAD, (sqrt (r2max) + 1.1)) + call apseti (py, PYNVER, nvertices) + } + + call sfree (sp) + return (nvertices) +end + + +# AP_YMKPOLY -- Mark the coordinates of a polygon on the display device. + +int procedure ap_ymkpoly (py, im, id, x, y, max_nvertices, idflush) + +pointer py # polyphot structure +pointer im # the input image descriptor +pointer id # display pointer +real x[ARB] # x coords of vertices +real y[ARB] # y coords of vertices +int max_nvertices # maximum number of vertices +int idflush # flush the imd interface ? + +int i, nvertices, stat, wcs, key +pointer sp, cmd +real xtemp, ytemp, xc, yc, r2max, r2 +real apstatr(), asumr() +int clgcur(), apstati() +errchk gscur + +begin + # Reopen the device. + if (id != NULL) + call greactivate (id, 0) + + # Initialize. + call smark (sp) + call salloc (cmd, SZ_LINE, TY_CHAR) + + # Type prompt. + call printf ( + "Mark polygon vertex [space=mark,q=quit].\n") + stat = clgcur ("icommands", xtemp, ytemp, wcs, key, Memc[cmd], SZ_LINE) + call ap_vtol (im, xtemp, ytemp, xtemp, ytemp, 1) + + + # Fetch the polygon and draw it on the display. + nvertices = 0 + while (stat != EOF) { + + if (key == 'q') + break + + # Decode and draw vertices. + nvertices = nvertices + 1 + if (nvertices <= max_nvertices) { + x[nvertices] = xtemp + y[nvertices] = ytemp + if (id != NULL) { + if (nvertices == 1) + call gamove (id, x[1], y[1]) + else { + call gadraw (id, x[nvertices], y[nvertices]) + if (idflush == YES) + call gframe (id) + else + call gflush (id) + } + } + } else + break + + # Type prompt. + call printf ( + "Mark polygon vertex [space=mark,q=quit].\n") + stat = clgcur ("icommands", xtemp, ytemp, wcs, key, Memc[cmd], + SZ_LINE) + call ap_vtol (im, xtemp, ytemp, xtemp, ytemp, 1) + } + call printf ("\n") + + call sfree (sp) + + # Return EOF or the number of vertices in the polygon. + if (stat == EOF) { + + call apsetr (py, PYXMEAN, INDEFR) + call apsetr (py, PYYMEAN, INDEFR) + call apsetr (py, OPYXMEAN, INDEFR) + call apsetr (py, OPYYMEAN, INDEFR) + call apsetr (py, PYCX, INDEFR) + call apsetr (py, PYCY, INDEFR) + call apsetr (py, OPYCX, INDEFR) + call apsetr (py, OPYCY, INDEFR) + call apsetr (py, PYMINRAD, INDEFR) + call apseti (py, PYNVER, 0) + nvertices = EOF + + } else if (nvertices <= 2) { + + call apsetr (py, PYXMEAN, INDEFR) + call apsetr (py, PYYMEAN, INDEFR) + call apsetr (py, OPYXMEAN, INDEFR) + call apsetr (py, OPYYMEAN, INDEFR) + call apsetr (py, PYCX, INDEFR) + call apsetr (py, PYCY, INDEFR) + call apsetr (py, OPYCX, INDEFR) + call apsetr (py, OPYCY, INDEFR) + call apsetr (py, PYMINRAD, INDEFR) + call apseti (py, PYNVER, 0) + nvertices = 0 + + } else { + + # Add the last vertex and close the polygon. + x[nvertices+1] = x[1] + y[nvertices+1] = y[1] + if (id != NULL) { + call gadraw (id, x[1], y[1]) + if (idflush == YES) + call gframe (id) + else + call gflush (id) + } + + # Compute and save the mean polygon coords. + xc = asumr (x, nvertices) / nvertices + yc = asumr (y, nvertices) / nvertices + call apsetr (py, PYXMEAN, xc) + call apsetr (py, PYYMEAN, yc) + call apsetr (py, PYCX, xc) + call apsetr (py, PYCY, yc) + if (id != NULL) + call gscur (id, apstatr (py, PYCX), apstatr (py, PYCY)) + + # Compute the mean sky annulus. + r2max = 0.0 + do i = 1, nvertices { + r2 = (x[i] - xc) ** 2 + (y[i] - yc) ** 2 + if (r2 > r2max) + r2max = r2 + } + + # Compute output coordinate centers. + switch (apstati(py,WCSOUT)) { + case WCS_WORLD, WCS_PHYSICAL: + call ap_ltoo (py, xc, yc, xc, yc, 1) + case WCS_TV: + call ap_ltov (im, xc, yc, xc, yc, 1) + default: + ; + } + call apsetr (py, OPYXMEAN, xc) + call apsetr (py, OPYYMEAN, yc) + call apsetr (py, OPYCX, xc) + call apsetr (py, OPYCY, yc) + + call apseti (py, PYNVER, nvertices) + call apsetr (py, PYMINRAD, (sqrt (r2max) + 1.1)) + + } + + if (id != NULL) + call gdeactivate (id, 0) + + return (nvertices) +end -- cgit