aboutsummaryrefslogtreecommitdiff
path: root/pkg/xtools/skywcs/skwrite.x
diff options
context:
space:
mode:
Diffstat (limited to 'pkg/xtools/skywcs/skwrite.x')
-rw-r--r--pkg/xtools/skywcs/skwrite.x510
1 files changed, 510 insertions, 0 deletions
diff --git a/pkg/xtools/skywcs/skwrite.x b/pkg/xtools/skywcs/skwrite.x
new file mode 100644
index 00000000..2e779b09
--- /dev/null
+++ b/pkg/xtools/skywcs/skwrite.x
@@ -0,0 +1,510 @@
+include "skywcsdef.h"
+include "skywcs.h"
+
+
+# SK_IIPRINT -- Print a summary of the input image or list coordinate system.
+
+procedure sk_iiprint (label, imagesys, mw, coo)
+
+char label[ARB] #I the input label
+char imagesys[ARB] #I the input image name and wcs
+pointer mw #I pointer to the image wcs
+pointer coo #I pointer to the coordinate system structure
+
+begin
+ if (mw == NULL)
+ call sk_inprint (label, imagesys, SKY_CTYPE(coo),
+ SKY_RADECSYS(coo), SKY_EQUINOX(coo), SKY_EPOCH(coo))
+ else
+ call sk_imprint (label, imagesys, SKY_CTYPE(coo), SKY_PLNGAX(coo),
+ SKY_PLATAX(coo), SKY_WTYPE(coo), SKY_PIXTYPE(coo),
+ SKY_RADECSYS(coo), SKY_EQUINOX(coo), SKY_EPOCH(coo))
+end
+
+
+# SK_IIWRITE -- Write a summary of the input image or list coordinate system
+# to the output file
+
+procedure sk_iiwrite (fd, label, imagesys, mw, coo)
+
+int fd #I the output file descriptor
+char label[ARB] #I the input label
+char imagesys[ARB] #I the input image name and wcs
+pointer mw #I pointer to the image wcs
+pointer coo #I pointer to the coordinate system structure
+
+begin
+ if (mw == NULL)
+ call sk_inwrite (fd, label, imagesys, SKY_CTYPE(coo),
+ SKY_RADECSYS(coo), SKY_EQUINOX(coo), SKY_EPOCH(coo))
+ else
+ call sk_imwrite (fd, label, imagesys, SKY_CTYPE(coo),
+ SKY_PLNGAX(coo), SKY_PLATAX(coo), SKY_WTYPE(coo),
+ SKY_PIXTYPE(coo), SKY_RADECSYS(coo), SKY_EQUINOX(coo),
+ SKY_EPOCH(coo))
+end
+
+
+# SK_INPRINT -- Print a summary of the input list coordinate system.
+# This should probably be a call to sk_inwrite with the file descriptor
+# set to STDOUT to avoid duplication of code. There was a reason for
+# having two routines at one point but I can't remember what it was ...
+
+procedure sk_inprint (label, system, ctype, radecsys, equinox, epoch)
+
+char label[ARB] #I the input label
+char system[ARB] #I the input system
+int ctype #I the input coordinate type
+int radecsys #I the input equatorial reference system
+double equinox #I the input equinox
+double epoch #I the input epoch of the observation
+
+pointer sp, radecstr
+double sl_epj(), sl_epb()
+int sk_wrdstr()
+
+begin
+ call smark (sp)
+ call salloc (radecstr, SZ_FNAME, TY_CHAR)
+
+ switch (ctype) {
+
+ case CTYPE_EQUATORIAL:
+ if (sk_wrdstr (radecsys, Memc[radecstr], SZ_FNAME,
+ EQTYPE_LIST) <= 0)
+ call strcpy ("FK5", Memc[radecstr], SZ_FNAME)
+ call strupr (Memc[radecstr])
+ call printf ("%s: %s Coordinates: equatorial %s\n")
+ call pargstr (label)
+ call pargstr (system)
+ call pargstr (Memc[radecstr])
+ switch (radecsys) {
+ case EQTYPE_GAPPT:
+ call printf (" MJD: %0.5f Epoch: J%0.8f B%0.8f\n")
+ call pargd (epoch)
+ if (IS_INDEFD(epoch)) {
+ call pargd (INDEFD)
+ call pargd (INDEFD)
+ } else {
+ call pargd (sl_epj (epoch))
+ call pargd (sl_epb (epoch))
+ }
+ case EQTYPE_FK5, EQTYPE_ICRS:
+ call printf (" Equinox: J%0.3f Epoch: J%0.8f MJD: %0.5f\n")
+ call pargd (equinox)
+ call pargd (sl_epj(epoch))
+ call pargd (epoch)
+ default:
+ call printf (" Equinox: B%0.3f Epoch: B%0.8f MJD: %0.5f\n")
+ call pargd (equinox)
+ call pargd (sl_epb(epoch))
+ call pargd (epoch)
+ }
+
+ case CTYPE_ECLIPTIC:
+ call printf ("%s: %s Coordinates: ecliptic\n")
+ call pargstr (label)
+ call pargstr (system)
+ call printf (" MJD: %0.5f Epoch: J%0.8f B%0.8f\n")
+ call pargd (epoch)
+ if (IS_INDEFD(epoch)) {
+ call pargd (INDEFD)
+ call pargd (INDEFD)
+ } else {
+ call pargd (sl_epj(epoch))
+ call pargd (sl_epb(epoch))
+ }
+
+ case CTYPE_GALACTIC:
+ call printf ("%s: %s Coordinates: galactic\n")
+ call pargstr (label)
+ call pargstr (system)
+ call printf (" MJD: %0.5f Epoch: J%0.8f B%0.8f\n")
+ call pargd (epoch)
+ call pargd (sl_epj (epoch))
+ call pargd (sl_epb (epoch))
+
+ case CTYPE_SUPERGALACTIC:
+ call printf ("%s: %s Coordinates: supergalactic\n")
+ call pargstr (label)
+ call pargstr (system)
+ call printf (" MJD: %0.5f Epoch: J%0.8f B%0.8f\n")
+ call pargd (epoch)
+ call pargd (sl_epj (epoch))
+ call pargd (sl_epb (epoch))
+
+ }
+
+ call sfree (sp)
+end
+
+
+# SK_INWRITE -- Write a summary of the input coordinate system.
+
+procedure sk_inwrite (fd, label, system, ctype, radecsys, equinox, epoch)
+
+int fd #I the output file descriptor
+char label[ARB] #I the input label
+char system[ARB] #I the input system
+int ctype #I the input coordinate type
+int radecsys #I the input equatorial reference system
+double equinox #I the input equinox
+double epoch #I the input epoch of the observation
+
+pointer sp, radecstr
+double sl_epj(), sl_epb()
+int sk_wrdstr()
+
+begin
+ call smark (sp)
+ call salloc (radecstr, SZ_FNAME, TY_CHAR)
+
+ switch (ctype) {
+
+ case CTYPE_EQUATORIAL:
+ if (sk_wrdstr (radecsys, Memc[radecstr], SZ_FNAME,
+ EQTYPE_LIST) <= 0)
+ call strcpy ("FK5", Memc[radecstr], SZ_FNAME)
+ call strupr (Memc[radecstr])
+ call fprintf (fd, "# %s: %s Coordinates: equatorial %s\n")
+ call pargstr (label)
+ call pargstr (system)
+ call pargstr (Memc[radecstr])
+ switch (radecsys) {
+ case EQTYPE_GAPPT:
+ call fprintf (fd, "# MJD: %0.5f Epoch: J%0.8f B%0.8f\n")
+ call pargd (epoch)
+ if (IS_INDEFD(epoch)) {
+ call pargd (INDEFD)
+ call pargd (INDEFD)
+ } else {
+ call pargd (sl_epj(epoch))
+ call pargd (sl_epb(epoch))
+ }
+ case EQTYPE_FK5, EQTYPE_ICRS:
+ call fprintf (fd,
+ "# Equinox: J%0.3f Epoch: J%0.8f MJD: %0.5f\n")
+ call pargd (equinox)
+ call pargd (sl_epj(epoch))
+ call pargd (epoch)
+ default:
+ call fprintf (fd,
+ "# Equinox: B%0.3f Epoch: B%0.8f MJD: %0.5f\n")
+ call pargd (equinox)
+ call pargd (sl_epb(epoch))
+ call pargd (epoch)
+ }
+
+ case CTYPE_ECLIPTIC:
+ call fprintf (fd, "# %s: %s Coordinates: ecliptic\n")
+ call pargstr (label)
+ call pargstr (system)
+ call fprintf (fd, "# MJD: %0.5f Epoch: J%0.8f B%0.8f\n")
+ call pargd (epoch)
+ if (IS_INDEFD(epoch)) {
+ call pargd (INDEFD)
+ call pargd (INDEFD)
+ } else {
+ call pargd (sl_epj(epoch))
+ call pargd (sl_epb(epoch))
+ }
+
+ case CTYPE_GALACTIC:
+ call fprintf (fd, "# %s: %s Coordinates: galactic\n")
+ call pargstr (label)
+ call pargstr (system)
+ call fprintf (fd, "# MJD: %0.5f Epoch: J%0.8f B%0.8f\n")
+ call pargd (epoch)
+ call pargd (sl_epj(epoch))
+ call pargd (sl_epb(epoch))
+
+ case CTYPE_SUPERGALACTIC:
+ call fprintf (fd, "# %s: %s Coordinates: supergalactic\n")
+ call pargstr (label)
+ call pargstr (system)
+ call fprintf (fd, "# MJD: %0.5f Epoch: J%0.8f B%0.8f\n")
+ call pargd (epoch)
+ call pargd (sl_epj(epoch))
+ call pargd (sl_epb(epoch))
+
+ }
+
+ call sfree (sp)
+end
+
+
+# SK_IMPRINT -- Print a summary of the input image coordinate system.
+# This should probably be a call to sk_imwrite with the file descriptor
+# set to STDOUT to avoid duplication of code. There was a reason for
+# having two routines at one point but I can't remember what it was ...
+
+procedure sk_imprint (label, imagesys, ctype, lngax, latax, wtype, ptype,
+ radecsys, equinox, epoch)
+
+char label[ARB] #I input label
+char imagesys[ARB] #I the input image name and system
+int ctype #I the image coordinate type
+int lngax #I the image ra/glon/elon axis
+int latax #I the image dec/glat/elat axis
+int wtype #I the image projection type
+int ptype #I the image image wcs type
+int radecsys #I the image equatorial reference system
+double equinox #I the image equinox
+double epoch #I the image epoch of the observation
+
+pointer sp, imname, projstr, wcsstr, radecstr
+double sl_epj(), sl_epb()
+int sk_wrdstr()
+
+begin
+ call smark (sp)
+ call salloc (imname, SZ_FNAME, TY_CHAR)
+ call salloc (projstr, SZ_FNAME, TY_CHAR)
+ call salloc (wcsstr, SZ_FNAME, TY_CHAR)
+ call salloc (radecstr, SZ_FNAME, TY_CHAR)
+
+ call sscan (imagesys)
+ call gargwrd (Memc[imname], SZ_FNAME)
+ if (sk_wrdstr (wtype, Memc[projstr], SZ_FNAME, WTYPE_LIST) <= 0)
+ call strcpy ("linear", Memc[projstr], SZ_FNAME)
+ call strupr (Memc[projstr])
+ if (sk_wrdstr (ptype, Memc[wcsstr], SZ_FNAME, PIXTYPE_LIST) <= 0)
+ call strcpy ("world", Memc[wcsstr], SZ_FNAME)
+ call strlwr (Memc[wcsstr])
+
+ switch (ctype) {
+
+ case CTYPE_EQUATORIAL:
+ if (sk_wrdstr (radecsys, Memc[radecstr], SZ_FNAME,
+ EQTYPE_LIST) <= 0)
+ call strcpy ("FK5", Memc[radecstr], SZ_FNAME)
+ call strupr (Memc[radecstr])
+ call printf (
+ "%s: %s %s Projection: %s Ra/Dec axes: %d/%d\n")
+ call pargstr (label)
+ call pargstr (Memc[imname])
+ call pargstr (Memc[wcsstr])
+ call pargstr (Memc[projstr])
+ call pargi (lngax)
+ call pargi (latax)
+ switch (radecsys) {
+ case EQTYPE_GAPPT:
+ call printf (" Coordinates: equatorial %s\n")
+ call pargstr (Memc[radecstr])
+ call printf (" MJD: %0.5f Epoch: J%0.8f B%0.8f\n")
+ call pargd (epoch)
+ if (IS_INDEFD(epoch)) {
+ call pargd (INDEFD)
+ call pargd (INDEFD)
+ } else {
+ call pargd (sl_epj(epoch))
+ call pargd (sl_epb(epoch))
+ }
+ case EQTYPE_FK5, EQTYPE_ICRS:
+ call printf (" Coordinates: equatorial %s Equinox: J%0.3f\n")
+ call pargstr (Memc[radecstr])
+ call pargd (equinox)
+ call printf (" Epoch: J%0.8f MJD: %0.5f\n")
+ call pargd (sl_epj (epoch))
+ call pargd (epoch)
+ default:
+ call printf (" Coordinates: equatorial %s Equinox: B%0.3f\n")
+ call pargstr (Memc[radecstr])
+ call pargd (equinox)
+ call printf (" Epoch: B%0.8f MJD: %0.5f\n")
+ call pargd (sl_epb (epoch))
+ call pargd (epoch)
+ }
+
+ case CTYPE_ECLIPTIC:
+ call printf (
+ "%s: %s %s Projection: %s Elong/Elat axes: %d/%d\n")
+ call pargstr (label)
+ call pargstr (Memc[imname])
+ call pargstr (Memc[wcsstr])
+ call pargstr (Memc[projstr])
+ call pargi (lngax)
+ call pargi (latax)
+ call printf (" Coordinates: ecliptic\n")
+ call printf (" MJD: %0.5f Epoch: J%0.8f B%0.8f\n")
+ call pargd (epoch)
+ if (IS_INDEFD(epoch)) {
+ call pargd (INDEFD)
+ call pargd (INDEFD)
+ } else {
+ call pargd (sl_epj(epoch))
+ call pargd (sl_epb(epoch))
+ }
+
+ case CTYPE_GALACTIC:
+ call printf (
+ "%s: %s %s Projection: %s Glong/Glat axes: %d/%d\n")
+ call pargstr (label)
+ call pargstr (Memc[imname])
+ call pargstr (Memc[wcsstr])
+ call pargstr (Memc[projstr])
+ call pargi (lngax)
+ call pargi (latax)
+ call printf (" Coordinates: galactic\n")
+ call printf (" MJD: %0.5f Epoch: J%0.8f B%0.8f\n")
+ call pargd (epoch)
+ call pargd (sl_epj (epoch))
+ call pargd (sl_epb (epoch))
+
+ case CTYPE_SUPERGALACTIC:
+ call printf (
+ "%s: %s %s Projection: %s Slong/Slat axes: %d/%d\n")
+ call pargstr (label)
+ call pargstr (Memc[imname])
+ call pargstr (Memc[wcsstr])
+ call pargstr (Memc[projstr])
+ call pargi (lngax)
+ call pargi (latax)
+ call printf (" Coordinates: supergalactic\n")
+ call printf (" MJD: %0.5f Epoch: J%0.8f B%0.8f\n")
+ call pargd (epoch)
+ call pargd (sl_epj (epoch))
+ call pargd (sl_epb (epoch))
+ }
+
+ call sfree (sp)
+end
+
+
+# SK_IMWRITE -- Write a summary of the image coordinate system to the
+# output file.
+
+procedure sk_imwrite (fd, label, imagesys, ctype, lngax, latax, wtype, ptype,
+ radecsys, equinox, epoch)
+
+int fd #I the output file descriptor
+char label[ARB] #I input label
+char imagesys[ARB] #I the input image name and wcs
+int ctype #I the image coordinate type
+int lngax #I the image ra/glon/elon axis
+int latax #I the image dec/glat/elat axis
+int wtype #I the image projection type
+int ptype #I the image image wcs type
+int radecsys #I the image equatorial reference system
+double equinox #I the image equinox
+double epoch #I the image epoch of the observation
+
+pointer sp, imname, projstr, wcsstr, radecstr
+double sl_epj(), sl_epb()
+int sk_wrdstr()
+
+begin
+ call smark (sp)
+ call salloc (imname, SZ_FNAME, TY_CHAR)
+ call salloc (projstr, SZ_FNAME, TY_CHAR)
+ call salloc (wcsstr, SZ_FNAME, TY_CHAR)
+ call salloc (radecstr, SZ_FNAME, TY_CHAR)
+
+ call sscan (imagesys)
+ call gargwrd (Memc[imname], SZ_FNAME)
+ if (sk_wrdstr (wtype, Memc[projstr], SZ_FNAME, WTYPE_LIST) <= 0)
+ call strcpy ("linear", Memc[projstr], SZ_FNAME)
+ call strupr (Memc[projstr])
+ if (sk_wrdstr (ptype, Memc[wcsstr], SZ_FNAME, PIXTYPE_LIST) <= 0)
+ call strcpy ("world", Memc[wcsstr], SZ_FNAME)
+ call strlwr (Memc[wcsstr])
+
+ switch (ctype) {
+
+ case CTYPE_EQUATORIAL:
+ if (sk_wrdstr (radecsys, Memc[radecstr], SZ_FNAME,
+ EQTYPE_LIST) <= 0)
+ call strcpy ("FK5", Memc[radecstr], SZ_FNAME)
+ call strupr (Memc[radecstr])
+ call fprintf (fd,
+ "# %s: %s %s Projection: %s Ra/Dec axes: %d/%d\n")
+ call pargstr (label)
+ call pargstr (Memc[imname])
+ call pargstr (Memc[wcsstr])
+ call pargstr (Memc[projstr])
+ call pargi (lngax)
+ call pargi (latax)
+ switch (radecsys) {
+ case EQTYPE_GAPPT:
+ call fprintf (fd, "# Coordinates: equatorial %s\n")
+ call pargstr (Memc[radecstr])
+ call fprintf (fd, "# MJD: %0.5f Epoch: J%0.8f B%0.8f\n")
+ call pargd (epoch)
+ if (IS_INDEFD(epoch)) {
+ call pargd (INDEFD)
+ call pargd (INDEFD)
+ } else {
+ call pargd (sl_epj(epoch))
+ call pargd (sl_epb(epoch))
+ }
+ case EQTYPE_FK5, EQTYPE_ICRS:
+ call fprintf (fd,
+ "# Coordinates: equatorial %s Equinox: J%0.3f\n")
+ call pargstr (Memc[radecstr])
+ call pargd (equinox)
+ call fprintf (fd, "# Epoch: J%0.8f MJD: %0.5f\n")
+ call pargd (sl_epj(epoch))
+ call pargd (epoch)
+ default:
+ call fprintf (fd,
+ "# Coordinates: equatorial %s Equinox: B%0.3f\n")
+ call pargstr (Memc[radecstr])
+ call pargd (equinox)
+ call fprintf (fd, "# Epoch: B%0.8f MJD: %0.5f\n")
+ call pargd (sl_epb (epoch))
+ call pargd (epoch)
+ }
+
+ case CTYPE_ECLIPTIC:
+ call fprintf (fd,
+ "# %s: %s %s Projection: %s Elong/Elat axes: %d/%d\n")
+ call pargstr (label)
+ call pargstr (Memc[imname])
+ call pargstr (Memc[wcsstr])
+ call pargstr (Memc[projstr])
+ call pargi (lngax)
+ call pargi (latax)
+ call fprintf (fd, "# Coordinates: ecliptic\n")
+ call fprintf (fd, "# MJD: %0.5f Epoch: J%0.8f B%0.8f\n")
+ call pargd (epoch)
+ if (IS_INDEFD(epoch)) {
+ call pargd (INDEFD)
+ call pargd (INDEFD)
+ } else {
+ call pargd (sl_epj(epoch))
+ call pargd (sl_epb(epoch))
+ }
+
+ case CTYPE_GALACTIC:
+ call fprintf (fd,
+ "# %s: %s %s Projection: %s Glong/Glat axes: %d/%d\n")
+ call pargstr (label)
+ call pargstr (Memc[imname])
+ call pargstr (Memc[wcsstr])
+ call pargstr (Memc[projstr])
+ call pargi (lngax)
+ call pargi (latax)
+ call fprintf (fd, "# Coordinates: galactic\n")
+ call fprintf (fd, "# MJD: %0.5f Epoch: J%0.8f B%0.8f\n")
+ call pargd (epoch)
+ call pargd (sl_epj(epoch))
+ call pargd (sl_epb(epoch))
+
+ case CTYPE_SUPERGALACTIC:
+ call fprintf (fd,
+ "# %s: %s %s Projection: %s Slong/Slat axes: %d/%d\n")
+ call pargstr (label)
+ call pargstr (Memc[imname])
+ call pargstr (Memc[wcsstr])
+ call pargstr (Memc[projstr])
+ call pargi (lngax)
+ call pargi (latax)
+ call fprintf (fd, "# Coordinates: supergalactic\n")
+ call fprintf (fd, "# MJD: %0.5f Epoch: J%0.8f B%0.8f\n")
+ call pargd (epoch)
+ call pargd (sl_epj(epoch))
+ call pargd (sl_epb(epoch))
+ }
+
+ call sfree (sp)
+end
generated by cgit (git 2.34.1) at 2025-09-20 08:25:05 -0400