aboutsummaryrefslogtreecommitdiff
path: root/noao/onedspec/smw/funits.x
diff options
context:
space:
mode:
Diffstat (limited to 'noao/onedspec/smw/funits.x')
-rw-r--r--noao/onedspec/smw/funits.x445
1 files changed, 445 insertions, 0 deletions
diff --git a/noao/onedspec/smw/funits.x b/noao/onedspec/smw/funits.x
new file mode 100644
index 00000000..e3c8ddf5
--- /dev/null
+++ b/noao/onedspec/smw/funits.x
@@ -0,0 +1,445 @@
+include <ctype.h>
+include <error.h>
+include <funits.h>
+
+
+# FUN_OPEN -- Open funits package
+# It is allowed to open an unknown funit type
+
+pointer procedure fun_open (funits)
+
+char funits[ARB] # Units string
+pointer fun # Units pointer returned
+
+begin
+ call calloc (fun, FUN_LEN, TY_STRUCT)
+ iferr (call fun_decode (fun, funits)) {
+ call fun_close (fun)
+ call erract (EA_ERROR)
+ }
+ return (fun)
+end
+
+
+# FUN_CLOSE -- Close funits package
+
+procedure fun_close (fun)
+
+pointer fun # Units pointer
+
+begin
+ call mfree (fun, TY_STRUCT)
+end
+
+
+# FUN_COPY -- Copy funits pointer
+
+procedure fun_copy (fun1, fun2)
+
+pointer fun1, fun2 # Units pointers
+
+begin
+ if (fun2 == NULL)
+ call malloc (fun2, FUN_LEN, TY_STRUCT)
+ call amovi (Memi[fun1], Memi[fun2], FUN_LEN)
+end
+
+
+# FUN_DECODE -- Decode funits string and set up funits structure.
+# The main work is done in FUN_DECODE1 so that the funits string may
+# be recursive; i.e. the funits string may contain other funits strings.
+
+procedure fun_decode (fun, funits)
+
+pointer fun # Units pointer
+char funits[ARB] # Units string
+
+bool streq()
+pointer sp, funits1, temp
+errchk fun_decode1, fun_ctranr
+
+begin
+ if (streq (funits, FUN_USER(fun)))
+ return
+
+ call smark (sp)
+ call salloc (funits1, SZ_LINE, TY_CHAR)
+ call salloc (temp, FUN_LEN, TY_STRUCT)
+
+ # Save a copy to restore in case of an error.
+ call fun_copy (fun, temp)
+
+ iferr (call fun_decode1 (fun, funits, Memc[funits1], SZ_LINE)) {
+ call fun_copy (temp, fun)
+ call sfree (sp)
+ call erract (EA_ERROR)
+ }
+
+ call sfree (sp)
+end
+
+
+# FUN_DECODE1 -- Decode funits string and set up funits structure.
+# Unknown funit strings are allowed.
+
+procedure fun_decode1 (fun, funits, funits1, sz_funits1)
+
+pointer fun # Units pointer
+char funits[ARB] # Units string
+char funits1[sz_funits1] # Secondary funits string to return
+int sz_funits1 # Size of secondary funits string
+
+int funmod, funtype
+int i, j, k, nscan(), strdic(), strlen()
+real funscale
+pointer sp, str
+
+int class[FUN_NUNITS]
+real scale[FUN_NUNITS]
+data class /FUN_FREQ,FUN_FREQ,FUN_FREQ,FUN_WAVE/
+data scale /FUN_J,FUN_FU,FUN_CGSH,FUN_CGSA/
+
+begin
+ call smark (sp)
+ call salloc (str, SZ_FNAME, TY_CHAR)
+
+ call strcpy (funits, Memc[str], SZ_FNAME)
+ call strlwr (Memc[str])
+ call sscan (Memc[str])
+ funtype = 0
+ funmod = 0
+ do i = 1, 2 {
+ call gargwrd (Memc[str], SZ_FNAME)
+ if (nscan() != i)
+ break
+
+ j = strdic (Memc[str], Memc[str], SZ_FNAME, FUN_DIC)
+ for (k=strlen(Memc[str]); k>0 &&
+ (IS_WHITE(Memc[str+k-1]) || Memc[str+k-1]=='\n'); k=k-1)
+ Memc[str+k-1] = EOS
+
+ if (j > FUN_NUNITS) {
+ if (funmod != 0)
+ break
+ funmod = j - FUN_NUNITS
+ } else {
+ funtype = j
+ break
+ }
+ }
+ i = nscan()
+ call gargr (funscale)
+ if (nscan() != i+1)
+ funscale = 1
+
+ if (funtype == 0) {
+ FUN_TYPE(fun) = 0
+ FUN_CLASS(fun) = FUN_UNKNOWN
+ FUN_LABEL(fun) = EOS
+ call strcpy (funits, FUN_UNITS(fun), SZ_UNITS)
+ } else {
+ FUN_TYPE(fun) = funtype
+ FUN_CLASS(fun) = class[funtype]
+ FUN_MOD(fun) = funmod
+ FUN_SCALE(fun) = scale[funtype] * funscale
+ FUN_LABEL(fun) = EOS
+ FUN_UNITS(fun) = EOS
+ call strcpy (funits, FUN_USER(fun), SZ_UNITS)
+ switch (funmod) {
+ case FUN_LOG:
+ call strcat ("Log ", FUN_LABEL(fun), SZ_UNITS)
+ case FUN_MAG:
+ call strcat ("Mag ", FUN_LABEL(fun), SZ_UNITS)
+ }
+ call strcat ("Flux", FUN_LABEL(fun), SZ_UNITS)
+ if (funscale != 1) {
+ call sprintf (FUN_UNITS(fun), SZ_UNITS, "%sx%.1g")
+ call pargstr (Memc[str])
+ call pargr (funscale)
+ } else {
+ call sprintf (FUN_UNITS(fun), SZ_UNITS, "%s")
+ call pargstr (Memc[str])
+ }
+ }
+
+ call sfree (sp)
+end
+
+
+# FUN_COMPARE -- Compare two funits
+
+bool procedure fun_compare (fun1, fun2)
+
+pointer fun1, fun2 # Units pointers to compare
+bool strne()
+
+begin
+ if (strne (FUN_UNITS(fun1), FUN_UNITS(fun2)))
+ return (false)
+ if (strne (FUN_LABEL(fun1), FUN_LABEL(fun2)))
+ return (false)
+ return (true)
+end
+
+
+# FUN_CTRANR -- Transform funits
+# Error is returned if the transform cannot be made
+
+procedure fun_ctranr (fun1, fun2, dun, dval, fval1, fval2, nvals)
+
+pointer fun1 # Input funits pointer
+pointer fun2 # Output funits pointer
+pointer dun # Input units pointer
+real dval[nvals] # Input dispersion values
+real fval1[nvals] # Input flux values
+real fval2[nvals] # Output flux values
+int nvals # Number of values
+
+int i
+real s, lambda
+pointer ang, un_open()
+bool fun_compare()
+errchk un_open, un_ctranr
+
+begin
+ if (fun_compare (fun1, fun2)) {
+ call amovr (fval1, fval2, nvals)
+ return
+ }
+
+ if (FUN_CLASS(fun1) == FUN_UNKNOWN || FUN_CLASS(fun2) == FUN_UNKNOWN)
+ call error (1, "Cannot convert between selected funits")
+
+ call amovr (fval1, fval2, nvals)
+
+ s = FUN_SCALE(fun1)
+ switch (FUN_MOD(fun1)) {
+ case FUN_LOG:
+ do i = 1, nvals
+ fval2[i] = 10. ** fval2[i]
+ case FUN_MAG:
+ do i = 1, nvals
+ fval2[i] = 10. ** (-0.4 * fval2[i])
+ }
+ switch (FUN_CLASS(fun1)) {
+ case FUN_FREQ:
+ do i = 1, nvals
+ fval2[i] = fval2[i] / s
+ case FUN_WAVE:
+ if (FUN_CLASS(fun2) != FUN_WAVE) {
+ s = s * FUN_VLIGHT
+ ang = un_open ("angstroms")
+ do i = 1, nvals {
+ call un_ctranr (dun, ang, dval[i], lambda, 1)
+ fval2[i] = fval2[i] / s * lambda**2
+ }
+ call un_close (ang)
+ } else {
+ do i = 1, nvals
+ fval2[i] = fval2[i] / s
+ }
+ }
+
+ s = FUN_SCALE(fun2)
+ switch (FUN_CLASS(fun2)) {
+ case FUN_FREQ:
+ do i = 1, nvals
+ fval2[i] = fval2[i] * s
+ case FUN_WAVE:
+ if (FUN_CLASS(fun1) != FUN_WAVE) {
+ s = s * FUN_VLIGHT
+ ang = un_open ("angstroms")
+ do i = 1, nvals {
+ call un_ctranr (dun, ang, dval[i], lambda, 1)
+ fval2[i] = fval2[i] * s / lambda**2
+ }
+ call un_close (ang)
+ } else {
+ do i = 1, nvals
+ fval2[i] = fval2[i] * s
+ }
+ }
+ switch (FUN_MOD(fun2)) {
+ case FUN_LOG:
+ do i = 1, nvals
+ fval2[i] = log10 (fval2[i])
+ case FUN_MAG:
+ do i = 1, nvals
+ fval2[i] = -2.5 * log10 (fval2[i])
+ }
+end
+
+
+# FUN_CHANGER -- Change funits
+# Error is returned if the conversion cannot be made
+
+procedure fun_changer (fun, funits, dun, dvals, fvals, nvals, update)
+
+pointer fun # Units pointer (may be changed)
+char funits[ARB] # Desired funits
+pointer dun # Dispersion units pointer
+real dvals[nvals] # Dispersion values
+real fvals[nvals] # Flux Values
+int nvals # Number of values
+int update # Update funits pointer?
+
+bool streq(), fun_compare()
+pointer fun1, fun_open()
+errchk fun_open, fun_ctranr
+
+begin
+
+ # Check for same funit string
+ if (streq (funits, FUN_USER(fun)))
+ return
+
+ # Check for error in funits string, or the same funits.
+ fun1 = fun_open (funits)
+ if (fun_compare (fun1, fun)) {
+ call strcpy (funits, FUN_USER(fun), SZ_UNITS)
+ call fun_close (fun1)
+ return
+ }
+
+ iferr {
+ call fun_ctranr (fun, fun1, dun, dvals, fvals, fvals, nvals)
+ if (update == YES)
+ call fun_copy (fun1, fun)
+ call fun_close(fun1)
+ } then {
+ call fun_close(fun1)
+ call erract (EA_ERROR)
+ }
+end
+
+
+# FUN_CTRAND -- Transform funits
+# Error is returned if the transform cannot be made
+
+procedure fun_ctrand (fun1, fun2, dun, dval, fval1, fval2, nvals)
+
+pointer fun1 # Input funits pointer
+pointer fun2 # Output funits pointer
+pointer dun # Input dispersion units pointer
+double dval[nvals] # Input dispersion values
+double fval1[nvals] # Input flux values
+double fval2[nvals] # Output flux values
+int nvals # Number of values
+
+int i
+double s, lambda
+pointer ang, un_open()
+bool fun_compare()
+errchk un_open, un_ctrand
+
+begin
+ if (fun_compare (fun1, fun2)) {
+ call amovd (fval1, fval2, nvals)
+ return
+ }
+
+ if (FUN_CLASS(fun1) == FUN_UNKNOWN || FUN_CLASS(fun2) == FUN_UNKNOWN)
+ call error (1, "Cannot convert between selected funits")
+
+ call amovd (fval1, fval2, nvals)
+
+ s = FUN_SCALE(fun1)
+ switch (FUN_MOD(fun1)) {
+ case FUN_LOG:
+ do i = 1, nvals
+ fval2[i] = 10. ** fval2[i]
+ case FUN_MAG:
+ do i = 1, nvals
+ fval2[i] = 10. ** (-0.4 * fval2[i])
+ }
+ switch (FUN_CLASS(fun1)) {
+ case FUN_FREQ:
+ do i = 1, nvals
+ fval2[i] = fval2[i] / s
+ case FUN_WAVE:
+ if (FUN_CLASS(fun2) != FUN_WAVE) {
+ s = s * FUN_VLIGHT
+ ang = un_open ("angstroms")
+ do i = 1, nvals {
+ call un_ctrand (dun, ang, dval[i], lambda, 1)
+ fval2[i] = fval2[i] / s * lambda**2
+ }
+ call un_close (ang)
+ } else {
+ do i = 1, nvals
+ fval2[i] = fval2[i] / s
+ }
+ }
+
+ s = FUN_SCALE(fun2)
+ switch (FUN_CLASS(fun2)) {
+ case FUN_FREQ:
+ do i = 1, nvals
+ fval2[i] = fval2[i] * s
+ case FUN_WAVE:
+ if (FUN_CLASS(fun1) != FUN_WAVE) {
+ s = s * FUN_VLIGHT
+ ang = un_open ("angstroms")
+ do i = 1, nvals {
+ call un_ctrand (dun, ang, dval[i], lambda, 1)
+ fval2[i] = fval2[i] * s / lambda**2
+ }
+ call un_close (ang)
+ } else {
+ do i = 1, nvals
+ fval2[i] = fval2[i] * s
+ }
+ }
+ switch (FUN_MOD(fun2)) {
+ case FUN_LOG:
+ do i = 1, nvals
+ fval2[i] = log10 (fval2[i])
+ case FUN_MAG:
+ do i = 1, nvals
+ fval2[i] = -2.5 * log10 (fval2[i])
+ }
+
+end
+
+
+# FUN_CHANGED -- Change funits
+# Error is returned if the conversion cannot be made
+
+procedure fun_changed (fun, funits, dun, dvals, fvals, nvals, update)
+
+pointer fun # Units pointer (may be changed)
+char funits[ARB] # Desired funits
+pointer dun # Input dispersion pointer
+double dvals[nvals] # Input dispersion values
+double fvals[nvals] # Flux values
+int nvals # Number of values
+int update # Update funits pointer?
+
+bool streq(), fun_compare()
+pointer fun1, fun_open()
+errchk fun_open, fun_ctrand
+
+begin
+
+ # Check for same funit string
+ if (streq (funits, FUN_USER(fun)))
+ return
+
+ # Check for error in funits string, or the same funits.
+ fun1 = fun_open (funits)
+ if (fun_compare (fun1, fun)) {
+ call strcpy (funits, FUN_USER(fun), SZ_UNITS)
+ call fun_close (fun1)
+ return
+ }
+
+ iferr {
+ call fun_ctrand (fun, fun1, dun, dvals, fvals, fvals, nvals)
+ if (update == YES)
+ call fun_copy (fun1, fun)
+ call fun_close(fun1)
+ } then {
+ call fun_close(fun1)
+ call erract (EA_ERROR)
+ }
+end
generated by cgit (git 2.34.1) at 2025-09-20 07:13:00 -0400