Документ взят из кэша поисковой машины. Адрес оригинального документа : http://www.stsci.edu/ftp/software/ASpec/documents/misc/sub_calls.asc Дата изменения: Mon Dec 2 20:07:10 1996 Дата индексирования: Sun Dec 23 12:23:23 2007 Кодировка:

AS_STRUCT (Mar96) source/as_obj AS_STRUCT (Mar96)



NAME
as_alloc - Allocate the Aspect object & initialize to defaults from the CL
as_free - Free the ASpec object
as_clinit - Initialize the ASpec object from CL parameters


DESCRIPTION


ROUTINE DETAILS

pointer procedure as_alloc()
RETURNS
pointer to Aspec object

procedure as_free (o)
ARGUMENTS
o (input: pointer)
The Aspec object pointer. NULL on return

procedure as_clinit(o)
ARGUMENTS
o (input: pointer)
The Aspec object pointer.



AS_ADD_RANGE (4Jan96) source/lib AS_ADD_RANGE (4Jan96)



NAME
as_add_range -- Add domain range to fitting parameters.


DESCRIPTION
Ranges are "exclusive". I.e. there can be no overlapping ranges.
Therefore, when a new range is specified, it any part of the range
overlaps with one already defined, it will become "merged" with the
existing range.


ROUTINE DETAILS

procedure as_add_range (o, low, high)


ARGUMENTS
o (input: pointer)
The ASpec object

low (input: double)
Low domain limits

high (input: double)
High domain limits



AS_DEL_RANGE_FROM_C (4Jan96) source AS_DEL_RANGE_FROM_C (4Jan96)



NAME
as_del_range_from_coord -- Delete a domain range


DESCRIPTION


ROUTINE DETAILS

procedure as_del_range_from_coor (o, x)


ARGUMENTS
o (input: pointer)
ASpec object

x (input: double)
Coordinate within domain to delete



AS_DEL_RANGE (Aug94) source AS_DEL_RANGE (Aug94)



NAME
as_del_range -- Delete a domain range


DESCRIPTION



ROUTINE DETAILS

procedure as_del_range (o, id)


ARGUMENTS
o (input: pointer)
ASpec object

id (input: int)
Range identification



AS_COMP (Jan95) source/component AS_COMP (Jan95)



NAME
cp_alloc -- Allocate a component structure.
cp_free -- Free a component structure.
cp_copy -- Make a copy of a spectral component & return a pointer.


DESCRIPTION


ROUTINE DETAILS

pointer procedure cp_alloc (cnam, ic)
ARGUMENTS
cnam[ARB] (input: char)
Component name

ic (input: int)
Running component index

RETURNS
cp (output: pointer)
Pointer to component object


procedure cp_free (cp)
ARGUMENTS
cp (input: pointer)
Pointer to component object


pointer procedure cp_copy (cmp)
ARGUMENTS
cmp (input: pointer)
Pointer to original component object

RETURNS
cp (output: pointer)
Pointer to new component object



CP_CV_UNITS (30Jan96) source/component CP_CV_UNITS (30Jan96)



NAME
cp_cv_units - Determine whether component unit conversion is necessary.
cp_tophot - Convert flux/dispersion attributes to PHOTLAM units.
cp_fromphot - Convert flux/dispersion values from PHOTLAM to specified units.


DESCRIPTION


ROUTINE DETAILS

procedure cp_cv_units (c, funit, dunit)
ARGUMENTS
c (input: pointer)
Component to convert

funit (input: int)
Flux units to convert to

dunit (input: int)
Dispersion units to convert to


procedure cp_tophot (cp, n_cp)
ARGUMENTS
cp[ARB] (input: pointer)
Array of components

n_cp (input: int)
Number of components


procedure cp_fromphot (cp, n_cp, funits, dunits)
ARGUMENTS
cp[ARB] (input: pointer)
Array of components

n_cp (input: int)
Number of components

funits, dunits (input: int)
Flux/ dispersion units for output



CP_DEBUG (30Jan96) source/component CP_DEBUG (30Jan96)



NAME
cp_examine - Routine to dereference elements of component structure.
cp_view - Routine to print/examine contents of component structure.


DESCRIPTION


ROUTINE DETAILS

procedure cp_examine (o)
ARGUMENTS
o (input: pointer)
Object pointer


procedure cp_view (cnam, type, contrib, index, version, active, npar, r_wave
, use_rw, dunit, funit, val, hi, lo, sig, stp, tol, cons)

ARGUMENTS
cnam[ARB] (input: char)
Component name

type (input: int)
Component type

contrib (input: int)
Flux contribution type

index (input: int)
Component ID

version (input: int)
Version number

active (input: bool)
Is component active for fit?

npar (input: int)
Size of arrays

r_wave (input: real)
Reference wavelength

use_rw (input: bool)
Use reference wavelength?

dunit, funit (input: int)
Dispersion, flux units

val[npar] (input: double)
Value array

hi[npar] (input: double)
upper limit

lo[npar] (input: double)
Lower limit

stp[npar] (input: real)
Step array

tol[npar] (input: real)
Tolerance array

sig[npar] (input: real)
Sigma array

cons[npar] (input: short)
Constraint array



CP_GET (Jan96) source/component CP_GET (Jan96)



NAME
cpg_param - Return parameter value from parameter name.
cpg_pindx - Return parameter index from parameter name.
cpg_npars - Get the number of parameters associated with the component type.
cpg_rwave - Fetch appropriate reference wavelength.
cpg_dict - Return parameter name dictionary from component type.


DESCRIPTION


ROUTINE DETAILS

double procedure cpg_param (cp, pname)
ARGUMENTS
cp (input: pointer)
Component object

pname[ARB] (input: char)
Parameter name


int procedure cpg_pindx (cp_type, pname)
ARGUMENTS
cp_type (input: int)
Component type

pname[ARB] (input: char)
Parameter name


int procedure cpg_npars (cp_type)
ARGUMENTS
cp_type (input: int)
Component type


procedure cpg_rwave (cp, ref_wave)
ARGUMENTS
cp (input: pointer)
Individual component object

ref_wave (output: double)
Reference wavelength


procedure cpg_dict (cp_type, dictionary, max_chars)
ARGUMENTS
cp_type (input: int)
Component type

dictionary[ARB] (output: char)
Parameter name dictionary

max_chars (input: int)
Size of dictionary string



CP_LIST (24Oct94) source CP_LIST (24Oct94)



NAME
cp_list -- Component list manipulation routines


DESCRIPTION
Along with 'as_comp.h', this set of routines define the semantics of
component lists. Below is describe the component list structure
along with the routines to manipulate a list.

A component list is implemented as an array, each element of which
contains a pointer to a single component. The routines provide a
method of creating and destroying the list, adding and removing
elements, and of finding an element.

Also implemented is the concept of the "current component". A
component is current if it has just been added or queried for.


MEMORY STRUCTURE
This section describes the memory structure used to store the
component list and any macros associated with access to the memory
structure. All variables/macros begin with 'CPL_'.

CPL_A_PTR(o)
This contains the pointer to the array of elements. See also
CPL_A below.

CPL_A(o,i)
This is the component whose index in the array is 'i'.

CPL_A_SZ(o)
This contains the size of the array. Since components can be
added/deleted dynamically, the size of the array is not
necessarily the same as the number of components stored in the
array.

CPL_N(o)
This holds the number of components currently stored in the
list.

CPL_CUR(o)
This contains the id of the "current" component. If there is
no current component, the value is INDEFI.

CPL_LAST_ID(o)
The identifier that was given to the last component added to
the list.

CPL_SZ
Size of the memory structure.

CPL_GROW
The amount the array will grow (if necessary) when components
are added.


ROUTINE OVERVIEW
The following routines are available for manipulating component
lists. Information contained in a component list may also be
accessed directly from the memory structure. However, for
adding/deleting elements from the list, it is strongly suggested
that the routines described below be utilized.

cpl_alloc -- Create a component list
cpl_free -- Destroy a component list
cpl_add -- Add a component to the list
cpl_del -- Delete a component from the list
cpl_copy -- Copy a component list
cpl_find -- Find the index of a component in the list, given its id
cpl_get_cmp -- Get a pointer to a component from the list, given its id
cpl_append -- Append a component list onto another one.
cpl_mod_id -- Modify the id's of the components in the list.


ROUTINE DETAILS
Below, each routine is described in detail.

pointer cpl_alloc (size)
Create a component list

ARGUMENTS
size (input: integer)
Initial number of components the list will hold. If
unknown, specify INDEFI.

RETURNS
A pointer to the component list structure. This pointer
should be used in all macros and routines which use the
component list.

cpl_free (o, free)
Destroy a component list.

ARGUMENTS
o (input/output: pointer)
The pointer to the component list structure to be
freed. On return, the value will be NULL.

free (input: bool)
If 'true', each component contained in the list will be
destroyed using the 'cp_free' routine. If false, the
components will not be destroyed. This is useful if
there are other structures which point to the same
components.

cpl_add (o, c, newid)
Add a component to the list.

ARGUMENTS
o (input: pointer)
The component list pointer.

c (input: pointer)
The pointer to the component to add to the list.

newid (input: bool)
If 'true', the component will be assinged a new
identifier. If false, the component will retain its
identifier.

cpl_del (o, id, free)
Remove a component from the list based on the component id.

ARGUMENTS
o (input: pointer)
The component list pointer.

id (input: int)
The id of the component to remove from the list. If
INDEFI, the current component will be removed. If
there is no current component, an error will be
generated.

free (input: bool)
If 'true', the component will be destroyed using the
'cp_free' routine. Else, only the list reference to
the component will be destroyed, leaving the component
itself intact. Useful if other memory structures refer
to the same component.

pointer cpl_copy (o, copy)
Create a new component list that is a copy of the given one.

ARGUMENTS
o (input: pointer)
The component list to be copied.

copy (input: bool)
If 'true', a copy of each component will be created
using the 'cp_copy' routine. The new list will then
refer to the copies of the components. Else, the new
list will contain the same pointers to the components.
Hence, the following scenario can occure: If a copy of
a list is made using 'copy = false', then the new list
is destroyed with 'free = true', the original list will
now have no components to point to.

RETURNS
A pointer to the new list is returned.

int cpl_find (o, id)
Get the index into the array of a component with the specified
id.

ARGUMENTS
o (input: pointer)
The component list to search.

id (input: int)
The id of the component to find. If INDEFI, the
current component will be returned. If there is no
current component, an error will be generated.

RETURNS
An index into the array pointing to the specified component.

pointer cpl_get_cmp (o, id)
Get the component from the list with the specified id.

ARGUMENTS
o (input: pointer)
The component list to search.

id (input: int)
The id of the component to retrieve. If INDEFI, the
current component will be returned. If there is no
current component, an error will be generated.

RETURNS
The pointer to the desired component.

cpl_append (s, d, copy)
Append a component list onto another component list.

ARGUMENTS
s (input: pointer)
The "source" component list to append onto the
"destination" component list.

d (input/output: pointer)
The "destination" component list which will have the
source component list appended to it on exit. If NULL,
a copy of the source component list is returned.

copy (input: bool)
If 'true', a copy of each component is made. If false,
only the pointers to the components are copied.

cpl_mod_id (o, compact, offset)
Modify the identifiers of the components and their constraints
in the list. There are two operations that can be done to
modify the components. First, the id's can be "compacted",
i.e. the id's are made consecutive, starting at 1. This is
selected by the 'compact' argument.

Secondly, an offset can be added to all id's. If compaction is
selected, the offset is applied after compaction occurs.

ARGUMENTS
o (input: pointer)
The component list to compact the id's for.

compact (input: bool)
'true' if the id's should be "compacted". This means
that the id's of all components in the list will be
such that the id's run from 1 to n where n is the
number of components. The constraint expressions of
the components will be modified accordingly to maintain
consistency.

offset (input: int)
An offset to add to all component id's. The constraint
expression will be modified accordingly to maintain
consistency. If compaction was also requested, the
id's are first compacted, then the offset applied.

int aindxi (v, a, n)
Find a value in an array and return its index position in the
array. If the value is not found, an error is generated.

ARGUMENTS
v (input: int)
The value to find in the array.

a (input: int[n])
The array to look for the value in.

n (input: int)
The number of elements in the array

RETURNS (int)
The index into the array where the value is located.



CP_POPULATE (Jan96) source/component CP_POPULATE (Jan96)



NAME
cp_populate -- Populate a spectral component with initial values


DESCRIPTION


ROUTINE DETAILS

procedure cp_populate (cp, funit, dunit, ref_wav, val, hi, lo, stp, tol
sig, cnstr, exprn, pnam, family)

ARGUMENTS
cp (input: pointer)
Pointer to component object

funit (input: int)
Flux units for input values

dunit (input: int)
Dispersion units for input values

ref_wav (input: real)
Reference wavelength

val[ARB] (input: double)
Parameter values

hi[ARB] (input: double)
Upper limits on values

lo[ARB] (input: double)
Lower limits on values

stp[ARB] (input: real)
Initial step size for fit

tol[ARB] (input: real)
Fractional tolerance for fit

sig[ARB] (input: real)
Sigma of fitted parameter

cnstr[ARB] (input: short)
Constraint flag

exprn[SZ_EXPRESN, ARB] (input: char)
Constraint expression

pnam[SZ_PNAME, ARB] (input: char)
Parameter names

family[ARB] (input: char)
Function family name



CP_SPEC (Feb96) source/component CP_SPEC (Feb96)



NAME
cp_sp_update - Retrieve a user-defined component from disk if needed.
cpg_spect - Fetch user spectrum: wavelength, flux, and size.
cpg_sprange - Fetch wavelength wavelength extremes of user spectrum.
cpg_spunits - Fetch dispersion/flux units of user spectrum.


DESCRIPTION


ROUTINE DETAILS

procedure cp_sp_update (cp)
ARGUMENTS
cp (input: pointer)
Component object


procedure cpg_spect (cp, indx, uwave, uflux, npts)
ARGUMENTS
cp (input: pointer)
Component object

indx (input: int)
Index of desired spectrum

uwave (output: pointer)
Pointer to wavelength array

uflux (output: pointer)
Pointer to flux array

npts (output: int)
Size of returned arrays


procedure cpg_sprange (cp, indx, wave_hi, wave_lo)
ARGUMENTS
cp (input: pointer)
Component object

indx (input: int)
Index of desired spectrum

wave_hi (output: double)
Smallest wavelength in array

wave_lo (output: double)
Largest wavelength in array


procedure cpg_spunits (cp, indx, dunit, funit)
ARGUMENTS
cp (input: pointer)
Component object

indx (input: int)
Index of desired spectrum

dunit (output: int)
Dispersion unit of desired spectrum

funit (output: int)
Flux unit of desired spectrum



CP_UFUNC (Feb96) source/component CP_UFUNC (Feb96)



NAME
cp_uf_init - Initialize user function file preferences
uf_set_prefs - Set user function file preferences


DESCRIPTION


ROUTINE DETAILS

procedure cp_uf_init (as, cpl)
ARGUMENTS
as (input: pointer)
Aspec object

cpl (input: pointer)
Component list


procedure uf_set_prefs (fp)
ARGUMENTS
fp (input: pointer)
Component family attributes parameter object

s (input: pointer)
Spectrum parameter object



CP_UTIL (Jan96) source/component CP_UTIL (Jan96)



NAME
cpl_set_constr - Set constraint flag for all parameters in a component list
cp_set_contrib - Set flux contribution type for a component
cp_setid - Set a component id
cpl_bounds - Ensure valid parameter upper/lower bounds for component list.
cp_bounds - Ensure valid component parameter upper/lower bounds.


DESCRIPTION


ROUTINE DETAILS

procedure cpl_set_constr (cpl, value)
ARGUMENTS
cpl (input: pointer)
component list

value (input: int)
value of flag


procedure cps_contrib (cp)
ARGUMENTS
cp (input: pointer)
component object


procedure cp_setid (c, id)
ARGUMENTS
c (input: pointer)
component to set Id for

id (input: int)
Id to set


procedure cpl_bounds (cl, active_only)
ARGUMENTS
cl (input: pointer)
component list

active_only (input: bool)
affect only active components in list?


procedure cp_vounds (cp)
ARGUMENTS
cp (input: pointer)
component object



CP_VALID (Jan96) source/component CP_VALID (Jan96)



NAME
cp_validate - Validate parameter types for an input component.
cp_par_valid - Return ERR if any parameter has a value of INDEF.


DESCRIPTION


ROUTINE DETAILS

procedure cp_validate (cindx, npar, p_name, nchars, pindx)
ARGUMENTS
cindx (input: int)
index of component type

npar (input: int)
number of parameters

p_name[nchars, ARB] (input: char)
parameter names

nchars (input: int)
size of character strings

pindx[ARB] (output: int)
index of parameter dictionary


procedure cp_par_valid (cp)
ARGUMENTS
cp (input :pointer)
component object

RETURNS
status (output: int)
return status



CPL_CV_UNITS (09Jan96) source/component CPL_CV_UNITS (09Jan96)



NAME
cpl_cv_units - Convert a list of components to the specified units.


DESCRIPTION


ROUTINE DETAILS

procedure cpl_cv_units (l, funit, dunit)
ARGUMENTS
l (input: pointer)
component list to convert

funit (input: int)
new flux units

dunit (input: int)
new dispersion units



AS_DB (Jan96) source/database AS_DB (Jan96)



NAME
db_alloc -- Allocate the ASpect DataBase structure.
db_free -- Close the DataBase table & free the memory structure.
dbg_hkeyvals -- Fetch header keyword values from DataBase table.
db_col_init -- Fetch ASpect DataBase table pointers.
db_col_def -- Initialize the ASpect DataBase column attributes.


DESCRIPTION


ROUTINE DETAILS



AS_DBC (Apr95) source/database AS_DBC (Apr95)



NAME
dbc_add -- Write a new component to the DataBase table.
dbc_delete -- Delete a component from the DataBase table.
dbc_update -- Update the values of a component in the DataBase
table. dbc_get -- Fetch a component from the DataBase table.



AS_DBH (Aug94) source/database AS_DBH (Aug94)



NAME
dbh_geti -- Fetch header value TY_INT from the ASpect DataBase
table. dbh_gett -- Fetch header text keyword from the ASpect
DataBase table. dbh_seti -- Set/add header information, TY_INT, in
the ASpect DataBase table. dbh_sett -- Set/add header text
information in the ASpect DataBase table.



DB_UTIL (Dec95) source/database DB_UTIL (Dec95)



NAME
db_set_units - Set flux & dispersion units for a database object.
db_get_ncmp - Get no. components associated with a database object.



DBIO (Apr95) source/database DBIO (Apr95)



NAME
db_read -- Fetch ASpect component entries from the open DB table.
db_write -- Write out solutions to the open DataBase table.
db_index -- Construct an index for all components in the DataBase
table



AS_CONSTR (Jan96) source/fitting AS_CONSTR (Jan96)



NAME
as_constr - Evaluate the contraint expressions for all components.
aev_cfunc - Evaluate the named user function aev_getop - Called by
evvexpr to fetch operand


DESCRIPTION


ROUTINE DETAILS

double procedure as_constr (expr, cl)
ARGUMENTS
expr[ARB] (input: char)
constraint esxpression

cl (input: pointer)


procedure aev_getop (ed, opname, o)
ARGUMENTS
ed (input: pointer)
pointer to expression descriptor

opname[ARB] (input: char)
operand name

o (output: pointer)
operand pointer


procedure aev_cfunc (ed, fcname, args, nargs, o)
ARGUMENTS
ed (input: ponter)
pointer to expression object

fcname[ARB] (input: char)
name of user function to be executed

args[ARB] (input: pointer)
pointers to function arguments

nargs (input: int)
no. function arguments

o (output: pointer)
operand pointer



AS_FIT (Jan96) source/fitting AS_FIT (Jan96)



NAME
fit_alloc - Allocate the FIT object.
fit_copy - Make a copy of a FIT object.
fit_free - Deallocate memory & FIT object.


DESCRIPTION


ROUTINE DETAILS

pointer procedure fit_alloc ()
ARGUMENTS

RETURNS
pointer to fit object


pointer procedure fit_copy (ft)
ARGUMENTS
ft (input: pointer)
fit object to copy

RETURNS
new copy of input fit object


procedure fit_free (ft)
ARGUMENTS
ft (input: pointer)
fit object



AS_SIGMAS (Feb96) source/fitting AS_SIGMAS (Feb96)



NAME
do_errors - Derive fitted parameter errors with the specified method


DESCRIPTION


ROUTINE DETAILS

procedure do_errors (ft, merit)
ARGUMENTS
ft (input: pointer)
pointer to FIT control structure

merit (input: extern)
merit function to be evaluated



CHISPEC (Jan96) source/fitting CHISPEC (Jan96)



NAME
chispec - Evaluate the chi**2 of all fitted components to multiple spectra.
update - Updates free & constrained parameter values in active component.
set_lim - Enforces hard upper/lower limits on active component parameter values.


DESCRIPTION


ROUTINE DETAILS

real procedure chispec (ft, val)
ARGUMENTS
ft (input: pointer)
fit object

val[ARB] (input: double)
values of freely varying parameters


procedure update (cl, val)
ARGUMENTS
cl (input: pointer)
component list

val[ARB] (input: double)
values of freely varying parameters


procedure set_lim (cl)
ARGUMENTS
cl (input: pointer)
component list



CURV_MATRIX (Feb96) source/fitting CURV_MATRIX (Feb96)



NAME
err_from_curv_matrix - Evaluate parameter errors using a curvature matrix
geterrmat - Get diagonal elements of error matrix
curvmat - Calculate the curvature matrix


DESCRIPTION


ROUTINE DETAILS

procedure geterrmat (ft, errmat, np, merit)
ARGUMENTS
ft (input: pointer)
pointer to FIT control structure

errmat (input/output: double)
error matrix

np (input: int)
no. free parameters

merit (input: extern)
merrit function to be evaluated


procedure curvmat (ft, np, alpha, beta, delta, merit)
ARGUMENTS
ft (input: pointer)
pointer to FIT control structure

np (input: int)
no. free parameters

alpha[np, ARB] (output: double)
curvature of matrix

beta[ARB] (output: double)

delta[ARB] (output: double)

merit (input: extern)
merit function to be evaluated



EVAL_FIT (Jan96) source/fitting EVAL_FIT (Jan96)



NAME
ft_eval_comp - Evaluate components over a specified dispersion
range.


DESCRIPTION


ROUTINE DETAILS

pointer procedure ft_eval_comp (cl, d1, d2, npts, disp, dvec)
ARGUMENTS
cl (input: pointer)
pointer to component list

d1 (input: double)
first dispersion coordinates for fit

d2 (input: double)
last dispersion coordinates for fit

npts (input: int)
number of points evaluated in spectrum

disp (input: double)
dispersion

dvec[npts] (input: double)
dispersion vector

RETURNS
pointer to spectrum object



FT_CP_INIT (Jan96) source/fitting FT_CP_INIT (Jan96)



NAME
cpf_select - Select components to use in fit, and return no. free parameters
is_cp_active - Determine whether a component should be active during a fit.


DESCRIPTION


ROUTINE DETAILS

int procedure cpf_select (cl, fit_cl, domain, nfi, n_free)
ARGUMENTS
cl (input: pointer)
component list)

fit_cl (output: pointer)
list of components for use in fit

domain (input: pointer)
dispersion domains (intervals) for fit

nfi (input: int)
no. fit intervals

n_free (output: int)
no. free parameters


bool procedure is_cp_active (cp, domain, nfi)
ARGUMENTS
cp (input: pointer)
component descriptor

domain[2, ARB] (input: double)
array of dispersion domains for fit

nfi (input: int)
np.intervals in fit



FT_DOMAINS (Jan96) source/fitting FT_DOMAINS (Jan96)



NAME
zero_domain - Reset the domain counter to zero.
next_interval - Fetch the next domain reference spectrum & start/stop pixel.


DESCRIPTION


ROUTINE DETAILS

procedure zero_domain (ft)
ARGUMENTS
ft (input: pointer)
fit object


int procedure next_interval (ft, i_spec, first_pix, last_pix)
ARGUMENTS
ft (input: pointer)
the fit object

i_spec (output: int)
the fit object

first_pix (output: int)
first spectrum reference pixel

last_pix (output: int)
last spectrum reference pixel



FT_EVAL (Jan96) source/fitting FT_EVAL (Jan96)



NAME
eval_model - Evaluate model over wavelength array for all components.
eval_comp - Evaluate an individual component at wavelength X.
merit_units - Convert units to PHOTLAM for merit fn evaluation if needed.


DESCRIPTION


ROUTINE DETAILS

procedure eval_model (cl, x, val, npts)
ARGUMENTS
cl (input: pointer)
component list

x[ARB] (input: double)
wave array over which to evaluate components

val[ARB] (output: real)
function value @X

npts (input: int)
size of evaluation array


procedure eval_comp (cp, cp_type, x, val, npts)
ARGUMENTS
cp (input: pointer)
component descriptor

cp_type (input: int)
component function type

x[ARB] (input: double)
array of evaluation wavelengths

val[ARB] (output: real)
evaluated component array

npts (input: int)
size of evaluation array


procedure merit_units (cp, x, val, npts)
ARGUMENTS
cp (input: pointer)
component descriptor

x[ARB] (input: double)
array of evaluation wavelengths

val[ARB] (input/output: real)
evaluated component array

npts (input: int)
size of evaluation array



FT_EXE (Jan96) source/fitting FT_EXE (Jan96)



NAME
do_fit - Call the appropriate fit routine.


DESCRIPTION


ROUTINE DETAILS

int procedure do_fit (ft, report)
ARGUMENTS
ft ( input: pointer)
pointer to fit control structure

report (input: extern)
routine to call to report iterations



FT_INIT (Jan96) source/fitting FT_INIT (Jan96)



NAME
fit_init - Populate the FIT control memory structure. fit_seti -
Set individual parameters in FIT structure, TY_INT. fit_setr -
Set individual parameters in FIT structure, TY_REAL.


DESCRIPTION


ROUTINE DETAILS

procedure fit_init (ft, sl, cl, domain, nf_intervals, fit_typ, n_iter, ftol)
ARGUMENTS
ft (input: pointer)
fit object

sl (input: pointer)
list of spectrum objects

cl (input: pointer)
list of components

domain (input: pointer)
dispersion domains (intervals) for fit

nf_intervals (input: int)
no. fit intervals

fit_typ (input: int)
choice of fit algorithm

n_iter (input: int)
iteration limit for fit

ftol (input: real)
fractional tolerance for fit


procedure fit_seti (ft, param, value)
ARGUMENTS
ft (input: pointer)
fit object

param (input: pointer)
parameter index

value (input: int)
parameter value


procedure fit_setr (ft, param, value)
ARGUMENTS
ft (input: pointer)
fit object

param (input: pointer)
parameter index

value (input: int)
parameter value



FT_LIST (16Jun95) source FT_LIST (16Jun95)



NAME
ft_list -- Fit object list manipulation routines


DESCRIPTION
Along with 'as_fit.h', this set of routines define the semantics of
fit object lists. Below is describe the fit list structure along
with the routines to manipulate a list.

A fit list is implemented as an array, each element of which
contains a pointer to a single fit. The routines provide a method
of creating and destroying the list, adding and removing elements,
and of finding an element.

Also implemented is the concept of the "current fit". A fit is
current if it has just been added or queried for.


MEMORY STRUCTURE
This section describes the memory structure used to store the fit
list and any macros associated with access to the memory
structure. All variables/macros begin with 'FTL_'.

FTL_A_PTR(o)
This contains the pointer to the array of elements. See also
FTL_A below.

FTL_A(o,i)
This is the fit whose index in the array is 'i'.

FTL_A_SZ(o)
This contains the size of the array. Since spectra can be
added/deleted dynamically, the size of the array is not
necessarily the same as the number of spectra stored in the
array.

FTL_N(o)
This holds the number of spectra currently stored in the list.

FTL_CUR(o)
This contains the index into the array of the "current" fit. If
there is no current fit, the value is INDEFI.

FTL_LAST_ID(o)
The last id assigned to a fit in the list.

FTL_SZ
Size of the memory structure.

FTL_GROW
The amount the array will grow (if necessary) when spectra are
added.


ROUTINE OVERVIEW
The following routines are available for manipulating fit lists.
Information contained in a fit list may also be accessed directly
from the memory structure. However, for adding/deleting elements
from the list, it is strongly suggested that the routines described
below be utilized.

ftl_alloc -- Create a fit list
ftl_free -- Destroy a fit list
ftl_add -- Add a fit to the list
ftl_del -- Delete a fit from the list
ftl_copy -- Copy a fit list
ftl_find -- Find the index of a fit in the list, given its id
ftl_get_fit -- Get a pointer to a fit from the list, given its id


ROUTINE DETAILS
Below, each routine is described in detail.

pointer ftl_alloc (size)
Create a fit list

ARGUMENTS
size (input: integer)
Initial number of spectra the list will hold. If
unknown, specify INDEFI.

RETURNS
A pointer to the fit list structure. This pointer should be
used in all macros and routines which use the fit list.

ftl_free (o, free)
Destroy a fit list.

ARGUMENTS
o (input/output: pointer)
The pointer to the fit list structure to be freed. On
return, the value will be NULL.

free (input: bool)
If TRUE, each fit contained in the list will be
destroyed using the 'spunmap' routine. If false, the
spectra will not be destroyed. This is useful if there
are other structures which point to the same spectra.

ftl_add (o, s, newid)
Add a fit to the list.

ARGUMENTS
o (input: pointer)
The fit list pointer.

s (input: pointer)
The pointer to the fit to add to the list.

newid (input: bool)
TRUE if a new id should be assigned to the component
being added.

ftl_del (o, id, free)
Remove a fit from the list based on the fit id.

ARGUMENTS
o (input: pointer)
The fit list pointer.

id (input: int)
The id of the fit to remove from the list. If INDEFI,
the current fit will be removed. If there is no
current fit, an error will be generated.

free (input: bool)
If TRUE, the fit will be destroyed using the 'fit_free'
routine. Else, only the list reference to the fit will
be destroyed, leaving the fit itself intact. Useful if
other memory structures refer to the same fit.

pointer ftl_copy (o, copy)
Create a new fit list that is a copy of the give one.

ARGUMENTS
o (input: pointer)
The fit list to be copied.

copy (input: bool)
If TRUE, a copy of each fit will be created using the
'fit_copy' routine. The new list will then refer to
the copies of the fit objects. Else, the new list will
contain the same pointers to the fits. Hence, the
following scenario can occure: If a copy of a list is
made using 'copy = false', then the new list is
destroyed with 'free = true', the original list will
now have no fits to which to point.

RETURNS
A pointer to the new list is returned.

int ftl_find (o, id)
Get the index into the array of a fit with the specified id.

ARGUMENTS
o (input: pointer)
The fit list to search.

id (input: int)
The id of the fit to find. If INDEFI, the current fit
will be returned. If there is no current fit, an error
will be generated.

RETURNS
An index into the array pointing to the specified fit.

pointer ftl_get_fit (o, id)
Get the fit from the list with the specified id.

ARGUMENTS
o (input: pointer)
The fit list to search.

id (input: int)
The id of the fit to retrieve. If INDEFI, the current
fit will be returned. If there is no current fit, an
error will be generated.

RETURNS
The pointer to the desired fit.



FT_SP_INIT (Feb96) source/fitting FT_SP_INIT (Feb96)



NAME
spf_select - Select spectra to use in fit, and return no. fit intervals
set_reg_limits - Set up fit region delimitors from domains
spf_err_valid - Populates error arrays for each spectrum in a list if undefined
w_near - Returns element K of array X closest to (but less than) X1


DESCRIPTION


ROUTINE DETAILS

int procedure spf_select( sl_in, sl_out, domain, nfi )
ARGUMENTS
sl_in (input: pointer)
spectrum list

sl_out (output: pointer)
list of selected spectra

domain[2, ARB] (input: double)
array of first, last dispersn coords for fit

nfi (input: int)
no. fit intervals


bool procedure spf_err_valid (sl)
ARGUMENTS
sl (input: pointer)
spectrum list


procedure set_reg_limits (ft, domain, nfi, nt_pix, nt_rej)
ARGUMENTS
ft (input: pointer)
fit object

domain[2, ARB] (input: double)
array of first, last dispersn coords for fit

nfi (input: int)
no. fit intervals

nt_pix (output: int)
n_pix in all relevant spectrum

nt_rej (output: int)
n_rejected_pix in all relevant spectrum


procedure spf_select (sl)
ARGUMENTS
sl (input: pointer)
spectrum list


int procedure w_near (x, x_refer, np)
ARGUMENTS
x[ARB] (input: double)
array

x_refer (input: double)
reference value

np (input; int)
no. elements in X



MARQUADT (Jan96) source/fitting MARQUADT (Jan96)



NAME
mrqt_init - Construct the memory structure for the Marquadt algorithm
mrqt_free - Deallocate memory & release the Marquadt memory structure
marquadt - Perform a fit to 1 or more components using the marquadt algorithm


DESCRIPTION


ROUTINE DETAILS

pointer procedure mrqt_init (ft)
ARGUMENTS
ft (input: pointer)
pointer to FIT control structure

RETURNS
pointer


procedure mrqt_free (ft, mx)
ARGUMENTS
ft (input: pointer)
pointer to FIT control structure

mx (input: pointer)
pointer to Marquadt data structure


int procedure marquadt (ft, mx, merit, report)
ARGUMENTS
ft (input: pointer)
pointer to FIT control structure

mx (input: pointer)
pointer to Marquadt data structure

merit (input: extern)
merit function to be evaluated

report (input: extern)
function to report running iterations



MATRIX (Jan96) source/fitting MATRIX (Jan96)



NAME
matinv - Invert a symmetric matrix and return its determinant


DESCRIPTION


ROUTINE DETAILS

double procedure matinv (a, np)
ARGUMENTS
a[np,ARB] (input/output: double)
matrix to invert in-place

np (input: int)
size of matrix



REDLAW (Jan96) nebular/lib REDLAW (Jan96)



NAME
gal_redlaw -- Calculates the ISEF for the Galaxy by Savage & Mathis
(1979) lmc_redlaw -- Calculates the ISEF for the LMC by Howarth
(1983) smc_redlaw -- Calculates the ISEF for the SMC by Prevot, et
al. (1984) ccm_redlaw -- Calculates the ISEF for the Galaxy by
Cardelli, Clayton & Mathis (1989)


DESCRIPTION
The following procedures return an evaluation of the interstellar
extinction function (ISEF) as published by various authors. Often
the ISEF is expressed as A(lambda)/E(B-V), and is normalized such
that the function is zero at the V passband, and 1.00 at B. These
functions are expressed in the literature as tabulated functions,
polynomial expressions, or both. The parameter is wavelength in
inverse microns, and the result is color excess in magnitudes. The
functions are defined or extended to about lambda=1000 Ang; the
form of the ISEF is not known shortward of that.

The routines here return A(lambda)/A(V) by renormaling the
extinction function to 0.0 at zero inverse wavelength, and 1.0 at
V. The ratio of total to selective absorption in the V passband,
R_v, is thought to average 3.10 for most lines of sight in the
Galaxy and in the LMC; it may be somewhat lower in the SMC. It is
the responsibility of the calling function to re-normalize (i.e.,
multiply) the returned value by the choice of R_v. Note that R_v
is not a parameter for any of the ISEF's here except for CCM.


ROUTINE DETAILS

procedure gal_redlaw (wave, extl, npts)
ARGUMENTS
wave[ARB] (input: double)
wavelength of interest

extl[ARB] (output: real)
extinction evaluation array

npts (input: int)
size of evaluation/wave arrays


procedure ccm_redlaw (wave, extl, npts, rv)
ARGUMENTS
wave[ARB] (input: double)
wavelength of interest

extl[ARB] (output: real)
extinction evaluation array

npts (input: int)
size of evaluation/wave arrays

rv (input: real)
ratio of total to selective extinction


procedure lmc_redlaw (wave, extl, npts)
ARGUMENTS
wave[ARB] (input: double)
wavelength of interest

extl[ARB] (output: real)
extinction evaluation array

npts (input: int)
size of evaluation/wave arrays


procedure smc_redlaw (wave, extl, npts)
ARGUMENTS
wave[ARB] (input: double)
wavelength of interest

extl[ARB] (output: real)
extinction evaluation array

npts (input: int)
size of evaluation/wave arrays



SIMPLEX (Jan96) source/fitting SIMPLEX (Jan96)



NAME
splx_init - Construct the memory structure for the simplex algorithm
splx_free - Deallocate memory and release simplex structure
simplex - Perform a fit to 1 or more components using the simplex algorithm
extremes - Determine indexes of highest/next-highest/lowest values in an array
vnew - Recompute the merit function for a trial simplex vertex


DESCRIPTION


ROUTINE DETAILS

pointer procedure splx_init (ft, merit)
ARGUMENTS
ft (input: pointer)
pointer to FIT control structure

merit (input: extern)
merit function to be evaluated

RETURNS
pointer


procedure splx_free (ft, sx)
ARGUMENTS
ft (input: pointer)
pointer to FIT control structure

sx (output: pointer)
pointer to simplex data structure


int procedure simplex (ft, sx, merit, report)
ARGUMENTS
ft (input: pointer)
pointer to FIT control structure

sx (output: pointer)
pointer to simplex data structure

merit (input: extern)
merit function to be evaluated

report (input: extern)
function to report running iterations


procedure extremes (x, np, hi, nhi, low)
ARGUMENTS
x[ARB] (input: real)
input array

np (input: int)
no. elements to be examined

hi (input: int)
index of highest values

nhi (input: int)
index of next-highest values

low (input: int)
index of lowest values


real procedure vnew (sx, ft, merit, hi, fac)
ARGUMENTS
ft (input: pointer)
pointer to FIT control structure

sx (output: pointer)
pointer to simplex data structure

merit (input: extern)
merit function to be evaluated

hi (input: int)
index of vertex w/worst residual

fac (input: real)
factor by which vertex is to be varied



DISP_UNIT (Mar96) source/lib DISP_UNIT (Mar96)



NAME
ang2any -- Convert dispersion values from Angstroms to specified units.
any2ang -- Convert dispersion values to Angstroms from specified units.
wav_freq -- Convert dispersion values: Angstroms to/from frequencies.
wav_ev -- Convert dispersion values: Angstroms to/from eV.


ROUTINE DETAILS

procedure ang2any (dunit, wave, dispers, npts)
ARGUMENTS
dunit (input: int)
Dispersion units on output.

wave[ARB] (input: double)
Wavelength array (Angstrom units)

dispers[ARB] (output: double)
Dispersion array (specified units)

npts (input: int)
Size of arrays.


procedure any2ang (dunit, dispers, wave, npts)
ARGUMENTS
dunit (input: int)
units of input dispersion

dispers[ARB] (input: double)
dispersion interval (frequency)

wave[ARB] (output: double)
wavelength in Angstroms

npts (input: int)
size of arrays


double procedure wav_freq (in)
ARGUMENTS
in (input: double)
Wavelength or frequency to be converted.

RETURNS
Converted frequency or wavelength.


double procedure wav_ev (in)
ARGUMENTS
in (input: double)
Wavelength or energy (in eV) to be converted.

RETURNS
Converted energy or wavelength.



DOMAINS (Mar96) source/lib DOMAINS (Mar96)



NAME
set_domain - Build a domain list from the fit_domains CL parameter.
next_domain - Fetch the next domain string from the domain list.
free_domain - Free the domain structure.
get_domains - Read a list of dispersion domains from a string.
incr_domain - Increment domain counter unless the dispersion domains overlap


DESRIPTION


ROUTINE DETAILS

pointer procedure set_domain (fdomain)
ARGUMENTS
fdomain[ARB] (input: char)
input fit_domains parameter value

RETURNS



int procedure next_domain (o, domain, max_chars)
ARGUMENTS
o (input: pointer)
the domain list object

domain[ARB] (output: char)
single domain string

max_chars (input: int)
size of domain string


procedure free_domain (o)
ARGUMENTS
o (input: pointer)
the domains object


int procedure get_domains (list, ds_domain, ds_min, ds_max, max_domains)
ARGUMENTS
list[ARB] (input: char)
string containing list of domains

ds_domain[2,ARB] (output: double)
array of dispersion domains

ds_min (input: double)
default dispersion minimum

ds_max (input: double)
default dispersion maximum

max_domains (input: int)
max permitted domains


int procedure incr_domain (ds_domain, indx)
ARGUMENTS
ds_domain[2,ARB] (input: double)
Array of dispersion domains

indx (input/output: int)
Index of current domain



GET_UNITS (Mar96) source/lib GET_UNITS (Mar96)



NAME
get_unit_typ - Interpret string to determine enumerated type of units.
get_unit_string - Get string corresponding to enumerated type of units.


DESCRIPTION


ROUTINE DETAILS

int procedure get_unit_typ (unit_string)
ARGUMENTS
unit_string[ARB] (input: char)
string containing units description


procedure get_unit_string (unit_type, unit_string, max_char)
ARGUMENTS
unit_type (input: int)
enumerated unit type

unit_string[ARB] (output: char)
string containing units description

max_char (input: int)
size of unit_string



INIT_CP (Jan96) source/lib INIT_CP (Jan96)



NAME
cp_init_from_marker - Generic component object initialization routine.
init_abedge - Initialize component object from marker WCS for ABEDGE
init_abgauss - Initialize component object from marker WCS for ABGAUSS
init_bpwrlaw - Initialize component object from marker WCS for BPWRLAW
init_exponen - Initialize component object from marker WCS for EXPONEN
init_gauss - Initialize component object from marker WCS for GAUSSIAN
init_linear - Initialize component object from marker WCS for LINEAR
init_lorentz - Initialize component object from marker WCS for LORENTZ
init_planck - Initialize component object from marker WCS for PLANCK
init_poly - Initialize component object from marker WCS for POLY
init_pwrlaw - Initialize component object from marker WCS for PWRLAW
init_voigt - Initialize component object from marker WCS for VOIGT


DESCRIPTION
These routines compute the initial values, upper/lower limits, etc.
of each parameter for a particular component type, given world
coordinate information that is derived from the GUI. It is the
responsibility of the calling routine to provide meaningful flux
and dispersion intervals, where d2, d1 and f2, f1 are the first and
second dispersion and flux values, respectively, in world
coordinates.


ROUTINE DETAILS

procedure cp_init_from_marker (cp, d1, d2, f1, f2, ftol)
ARGUMENTS
cp (input: pointer)
component descriptor

d1 (input: double)
minimum of dispersion interval in world coords

d2 (input: double)
maximum of dispersion interval in worl coords

f1 (input: real)
minimum of flux interval in world coordinates

f2 (input: real)
maximum of flux interval in world coordinates

ftol (input: real)
fractional fit tolerance for parameters


procedure init_abedge (cp, bx)
ARGUMENTS
cp (input: pointer)
component descriptor

bx (input: pointer)
box structure


init_abgauss (cp, bx)
ARGUMENTS
cp (input: pointer)
component descriptor

bx (input: pointer)
box structure


procedure init_bpwrlaw (cp, bx)
ARGUMENTS
cp (input: pointer)
component descriptor

bx (input: pointer)
box structure


procedure init_exponen (cp, bx)
ARGUMENTS
cp (input: pointer)
component descriptor

bx (input: pointer)
box structure


procedure init_gauss (cp, bx)
ARGUMENTS
cp (input: pointer)
component descriptor

bx (input: pointer)
box structure


procedure init_linear (cp, bx)
ARGUMENTS
cp (input: pointer)
component descriptor

bx (input: pointer)
box structure


procedure init_lorentz (cp, bx)
ARGUMENTS
cp (input: pointer)
component descriptor

bx (input: pointer)
box structure


procedure init_planck (cp, bx)
ARGUMENTS
cp (input: pointer)
component descriptor

bx (input: pointer)
box structure


procedure init_poly (cp, bx)
ARGUMENTS
cp (input: pointer)
component descriptor

bx (input: pointer)
box structure


procedure init_pwrlaw (cp, bx)
ARGUMENTS
cp (input: pointer)
component descriptor

bx (input: pointer)
box structure


procedure init_voigt (cp, bx)
ARGUMENTS
cp (input: pointer)
component descriptor

bx (input: pointer)
box structure



PARAMS (Mar96) source/lib PARAMS (Mar96)



NAME
params -- Parameter/value handling


DESCRIPTION
The following library handles the setting of arbitrary "parameters"
or variables.

There are two forms of input/output. The first is the direct
setting of a parameter or retrieval of a parameter. The second is
the use of strings for setting/inquiring. Specifically, a string
with the format:

param value param value ...

or alternatively with the restricted format:

param=value param=value ...

can be used to set parameter values, or to retrieve parameter
values. Note that parameter/value pairs are either separated by
spaces or by equals signs, depending upon how the "restrict_syntax"
value is set in the call to pr_alloc. If a string value contains a
space, the value should be quoted using double quotes. To a place
double quote within a value, escape it by preceeding the character
with a '\'. A backslash can be included by specifying "\\".

Note that all parameter values are stored as strings. Conversion
to/from strings are performed for all the type-dependent calls. The
parameter names are case-insensitive.


OBJECT MEMORY STRUCTURE
The following is the memory structure of the parameter object.
The primary variabls are listed first.

PNAMES
The string containing the list of parameters. The position
of a parameter name in the list determines where the value
of that parameter is in the VALUES array.

VALUE
The value of a parameter.

NPARS
The number of parameters defined.

The following variables/macros help manage the memory object.

PN_PTR
Pointer to the PNAMES string.

VA_PTR
Pointer to the V_PTR array.

V_PTR
The array of pointer to the values.

PN_SZ
Current maximum length of the PNAMES string

V_SZ
Current maximum number of values that can be stored.

PN_MAX_SZ
Current maximum length of a parameter name.

GROW
The amount the V_PTR array grows.

CURRENT
Index of the "current" parameter. This is the index of the
parameter returned by the last pr_nextX call.

RESTRICT_SYNTAX
Is the syntax of the param/value pairs to be restricted?
If so, the pairs must be separated with an equals sign.

OBJECT ROUTINE SUMMARY
The set of routines that make up the parameter object are
as follows:

pr_alloc -- Create a parameter object.
pr_free -- Free a parameter object.

pr_set -- Set parameter values from a parameter string.
pr_get -- Return parameter values in a string.
pr_getp -- Return parameter values as a string pointer.

pr_set[idsb] -- Set a specific parameter.
pr_get[idsb] -- Return value of a specific parameter.
pr_gsp -- Return pointer to value string.

pr_rew -- Rewind to beginning of list.
pr_next[idsb] -- Return value of the "next" parameter.
pr_nextsp -- Return pointer to value string of "next" parameter.

pr_copy -- Copy a params object.



STRCIC (Mar96) source/lib STRCIC (Mar96)



NAME
strcic -- Case-insensitive dictionary search


DESCRIPTION
Similar to strdic, except the string matching is case insensitive.
Both the dictionary and string to find are converted to lower case.
The string returned is how the string is found in the unmodified
dictionary.


ROUTINE DETAILS

int procedure strcic (in_str, out_str, maxchars, dict)
ARGUMENTS
in_str[ARB] (input: char)
input string to find in dictionary

out_str[maxchars] (output: char)
string found in dictionary

maxchars (input: int)
maximum size of out_str

dict[ARB] (input: char)
dictionary







NAME
str_to_d - Convert string to a number.


USAGE
nchar = ctod (str, ip, dval)


DESCRIPTION
Attempt to convert a string to a double precision number. This
routine differs from sys$fmtio/ctod.x in that it does not support
sexigesimal numbers. The index IP must be set to the first
character to be scanned upon entry to STR_TO_D, and will be left
pointing at the first untranslated character.

If the string is successfully converted, the number of characters
used is returned as the function argument. If the string (or the
first few characters of the string) cannot be interpreted as a
number, zero will be returned. Note that even if no numeric
characters are encountered, the index IP may be incremented, if
leading whitespace is encountered (but the return value N will
still be zero).

The upper case string "INDEF" is a legal real number, as is "." (.
== 0.0). Sexagesimal numbers are permitted. Excess digits of
precision are ignored. Out of range exponents are detected, and
result in the value INDEF being returned (this is not considered an
ERROR condition). Any number with an exponent greater than or
equal to MAX_EXPONENT is interpreted as INDEF, regardless of the
mantissa. The number need not contain a decimal point.

Sexagesimal numbers are not supported here because they are not
needed for this application, and because they interfere with using
a colon as a range separator.


ROUTINE DETAILS

int procedure str_to_d (str, ip, dval)
ARGUMENTS
str[ARB] (input: char)
string to be converted

ip (input/output: int)
pointer to a string

dval (output: double)
receives binary value



WORD (Mar96) source/lib WORD (Mar96)

NAME
asw_count -- Return the number of words in a list of words
asw_fetch -- Retrieve next word from string
asw_find -- Find the i-th word in a list of words
asw_match -- Return number of the word in the list which matches the word


DESCRIPTION
The procedures in this file perform simple processing on lists of
words. These procedures count the number of words in a list, fetch
the next word in a list, find the n-th word in a list, check for an
exact match between a word and a list of words. A word is any group
of contiguous characters which are neither whitespace or commas.
The definition of whitespace is anomalous, it includes any
character whose integer value is less than or equal to a blank.
Note that words cannot be delimeted by quotes and that escape
processing is not done.


ROUTINE DETAILS

int procedure asw_count (list)
ARGUMENTS
list[ARB] (input: char)
list of words


int procedure asw_fetch (str, ic, word, maxch)
ARGUMENTS
str[ARB] (input: char)
string containing words

ic (input/output: int)
index of starting character

word[ARB] (output: char)
word string

maxch (input: int)
declared length of output string


int procedure asw_find (index, list, word, maxch)
ARGUMENTS
index (input: int)
index to word within list

list[ARB] (input: char)
list of words

word[ARB] (output: char)
word returned by this procedure

maxch (input: int)
declared length of output string


int procedure asw_match (word, list)
ARGUMENTS
word[ARB] (input: char)
word to be matched

list[ARB] (input: char)
list of words



DEFCOL (20Dec95) source DEFCOL (20Dec95)



NAME
defcol -- Allocate the output spectrum columns


DESCRIPTION
Set up default format table for saving a spectrum.


ROUTINE DETAILS

procedure defcol (tp, otcp, colname, colunits, colfmt, cdtype, lendata, ncols)
ARGUMENTS
tp (input: pointer)
Table descriptor

otcp[ncols] (input: pointer)
Output column descriptors

colname[SZ_COLNAME,ncols] (input: char)
Column names

colunits[SZ_COLUNITS,ncols] (input: char)
Column units

colfmt[SZ_COLFMT,ncols] (input: int)
Column formats

cdtype[ncols] (input: int)
Column data types

lendata[ncols] (input: int)
Character sizes

ncols (input: int)
Number of columns to allocate



GET_COL_NAMES (15May96) source GET_COL_NAMES (15May96)



NAME
get_col_names -- Get spectrum table column names and default units
from the parameter object


DESCRIPTION
To extract the names of the columns from the parameter object, pass
the pointer of the parameter object to this routine and get back an
array of column names.


ROUTINE DETAILS

procedure get_col_names(pp,col_nam,col_unt)


ARGUMENTS
pp (input: pointer)
pointer to parameter object

col_nam (output: char)
array of column names

col_unt (output: char)
array of column units



GET_FILE_EXTS (15may96) source GET_FILE_EXTS (15may96)



NAME
get_file_exts -- Get input file extensions from the parameter object


DESCRIPTION
Pass the routine the pointer to the parameter object and scan for
the input file extensions.


ROUTINE DETAILS

procedure get_file_exts (pp,file_exts)


ARGUMENTS
pp (input: pointer)
pointer to parameter object

file_exts (output: char)
array of file extensions



GETHED (22Dec95) source GETHED (22Dec95)



NAME
gethead -- Get header parameters from input spectrum table


DESCRIPTION
This routines reads the values of the standard header keywords from
an input table. Supportted keywords are: TARGET, EXPTIME, TITLE,
GROUP, FLXIMAGE.


ROUTINE DETAILS

procedure gethed (tp, sp)


ARGUMENTS
tp (input: pointer)
Table descriptor

sp (input: pointer)
Spectrum object descriptor



GETIMD (22Dec95) GETIMD (22Dec95)



NAME
getimd -- Get header parameters from input spectrum image


DESCRIPTION
This routines reads the values of the standard header keywords from
an input image. Supported keywords are: TARGET, EXPTIME, TITLE,
GROUP, FLXIMAGE.



ROUTINE DETAILS

procedure get_im_head (im, sp)
ARGUMENTS


im (input: pointer)
Input image header

sp (input: pointer)
Input spectrum object



PUTHEAD (22Dec95) source PUTHEAD (22Dec95)



NAME
puthead -- Put header parameters in output spectrum table


DESCRIPTION
Add standard header parameters to an output table. Supported
keywords are: TARGET, EXPTIME, TITLE, GROUP, FLXIMAGE.


ROUTINE DETAILS

procedure puthead (sp, tp)


ARGUMENTS
sp (input: pointer)
Spectrum object descriptor

tp (input: pointer)
Table descriptor



SP_DEBUG (4Jan96) source SP_DEBUG (4Jan96)



NAME
sp_examine -- Routine to dereference arrays in spectrum structure for
examination w/the debugger.
sp_view -- Routine to print/examine contents of spectrum structure


DESCRIPTION
Routines to assist in debugging.


ROUTINE DETAILS

procedure sp_examine (o)
Dereference arrays in spectrum structure for examination w/the
debugger.

ARGUMENTS
o (input: pointer)
Object pointer

procedure sp_view (id, group, exptime, z_red, funit, dunit, np, d1, d2,
vers, f, w, s, m, sc, targ, title, fname)

Examine contents of spectrum structure.

ARGUMENTS
id (input: int)
Spectrum ID

group (input: int)
group of parent file

exptime (input: real)
Exposure time

z_red (input: real)
Redshift

funit (input: int)
Units of flux

dunit (input: int)
Dispersion arrays

np (input: int)
Size of arrays

d1,d2 (input: double)
First, last dispersion value

vers (input: int)
Spectrum version no.

f[np] (input: real)
Flux array

w[np] (input: double)
Dispersion array

s[np] (input: real)
Error array

m[np] (input: short)
Mask array

sc[np] (input: real)
Scale array

targ[ARB] (input: char)
Name of observed object

title[ARB] (input: char)
Title of observed object

fname[ARB] (input: char)
Name of parent file



SP_ALLOC (22Dec95) source SP_ALLOC (22Dec95)



NAME
sp_alloc -- Allocate space for spectrum object


DESCRIPTION
This routine creates an initialized spectrum object.


ROUTINE DETAILS

pointer procedure sp_alloc (npts)
ARGUMENTS
npts (input: int)
Number of data points in spectrum

RETURNS
sp (output: pointer)
Pointer to spectrum object



SP_COPY (1Feb96) source SP_COPY (1Feb96)



NAME
sp_copy -- Create a copy of a spectrum


DESCRIPTION
Create a duplicate of a spectrum object.


ROUTINE DETAILS

pointer procedure sp_copy (s)


ARGUMENTS
s (input: pointer)
Input spectrum to copy

RETURNS
n (output: pointer)
Pointer to new spectrum



SP_CV_UNITS (22Dec95) source SP_CV_UNITS (22Dec95)



NAME
sp_cv_units -- Convert spectrum to specified units


DESCRIPTION
Convert flux and wavelength vectors in a given spectrum objects
from current units to newly specified units for display purposes.


ROUTINE DETAILS

procedure sp_cv_units (s, funit, dunit)


ARGUMENTS
s (input: pointer)
Spectrum to convert

funit (input: int)
New flux units

dunit (input: int)
New dispersion units



SP_HST_OPEN (1Feb96) source SP_HST_OPEN (1Feb96)



NAME
sp_hst_open -- Fill a spectrum object from an HST format dataset


DESCRIPTION
Given a set of multigroup images with the same rootname and
different extensions(i.e., from GHRS or FOS), read in the spectrum
flux, wavelength, error, and mask data and fill a spectrum object.
The routine reads all available spectra and returns a pointer to
the list of newly read spectra.


ROUTINE DETAILS

pointer procedure sp_hst_open (ip, im, im_sec)


ARGUMENTS
ip (input: pointer)
pointer to CIF object

im (input: pointer)
pointer to primary CIF filename==flux image

im_sec[4] (input: pointer)
pointers to secondary CIF filenames

RETURNS
A pointer to a spectrum list of newly read spectra



SPIMAN (15May96) source SPIMAN (15May96)



NAME
sp_image_open -- Read an image (not a table) and fill a spectrum
object


DESCRIPTION
This routine reads a spectrum (or series of spectra) from a IRAF
image. Currently we can read:
OIF images with WCS information
STF multigroup images


ROUTINE DETAILS

pointer procedure sp_image_open (filename, file_exts)


ARGUMENTS
filename (input: char)
input filenam

file_exts (input: char)
default filename extensions

RETURNS
An integer value identifying the type of image.



SPIMAN (15May96) source SPIMAN (15May96)



NAME
sp_image_open -- Read an image (not a table) and fill a spectrum
object


DESCRIPTION
This routine reads a spectrum (or series of spectra) from a IRAF
image. Currently we can read:
OIF images with WCS information
STF multigroup images


ROUTINE DETAILS

pointer procedure sp_image_open (filename, file_exts)


ARGUMENTS
filename (input: char)
input filename

file_exts (input: char)
default filename extensions

RETURNS
A pointer to a list of spectra just read in



SPL_CV_UNITS (27Dec95) source SPL_CV_UNITS (27Dec95)



NAME
spl_cv_units -- Convert a list of spectra to the specified units


DESCRIPTION
This routine processes a list of spectra with successive calls to
sp_cv_units. For each spectrum in the list, the flux and wavelength
vectors are converted to the specified units.


ROUTINE DETAILS

procedure spl_cv_units (l, funit, dunit)


ARGUMENTS
l (input: pointer)
Pointer to spectrum list to convert

funit (input: int)
New flux units

dunit (input: int)
New dispersion units



SPL_DISP_DOMAIN (3Jan96) source/specio SPL_DISP_DOMAIN (3Jan96)



NAME
spl_disp_domain -- Find the dispersion extrema for the list of
spectra.


DESCRIPTION


ROUTINE DETAILS

procedure spl_disp_domain (s, low, high)


ARGUMENTS
s (input: pointer)
Spectrum list

low (output: double)
Low-valued dispersion

high (output: double)
High-valued dispersion



SPL_FLUX_DOMAIN (27Dec95) source SPL_FLUX_DOMAIN (27Dec95)



NAME
spl_flux_domain - Find flux range between two wavelengths for a spectrum list
sp_frange - Find flux range between two wavelengths for a single spectrum


DESCRIPTION


ROUTINE DETAILS

procedure spl_flux_domain (s, wlow, whigh, flow, fhigh)
ARGUMENTS
s (input: pointer)
Spectrum list

wlow (input: double)
Low dispersion

whigh (input: double)
High dispersion

flow (output: real)
Smallest flux

fhigh (output: real)
Largest flux

procedure sp_frange (s, wlow, whigh, flow, fhigh)
ARGUMENTS
spc (input: pointer)
Spectrum object

wlow (input: double)
Low dispersion

whigh (input: double)
High dispersion

flow (output: real)
Smallest flux

fhigh (output: real)
Largest flux



SP_LIST (16Jun95) source SP_LIST (16Jun95)



NAME
sp_list -- Spectrum list manipulation routines


DESCRIPTION
Along with 'as_spec.h', this set of routines define the semantics of
spectrum lists. Below is describe the spectrum list structure along
with the routines to manipulate a list.

A spectrum list is implemented as an array, each element of which
contains a pointer to a single spectrum. The routines provide a
method of creating and destroying the list, adding and removing
elements, and of finding an element.

Also implemented is the concept of the "current spectrum". A
spectrum is current if it has just been added or queried for.


MEMORY STRUCTURE
This section describes the memory structure used to store the
spectrum list and any macros associated with access to the memory
structure. All variables/macros begin with 'SPL_'.

SPL_A_PTR(o)
This contains the pointer to the array of elements. See also
SPL_A below.

SPL_A(o,i)
This is the spectrum whose index in the array is 'i'.

SPL_A_SZ(o)
This contains the size of the array. Since spectra can be
added/deleted dynamically, the size of the array is not
necessarily the same as the number of spectra stored in the
array.

SPL_N(o)
This holds the number of spectra currently stored in the list.

SPL_CUR(o)
This contains the index into the array of the "current"
spectrum. If there is no current spectrum, the value is
INDEFI.

SPL_LAST_ID(o)
The last id assigned to a spectrum in the list.

SPL_SZ
Size of the memory structure.

SPL_GROW
The amount the array will grow (if necessary) when spectra are
added.


ROUTINE OVERVIEW
The following routines are available for manipulating spectrum
lists. Information contained in a spectrum list may also be
accessed directly from the memory structure. However, for
adding/deleting elements from the list, it is strongly suggested
that the routines described below be utilized.

spl_alloc -- Create a spectrum list
spl_free -- Destroy a spectrum list
spl_add -- Add a spectrum to the list
spl_del -- Delete a spectrum from the list
spl_copy -- Copy a spectrum list
spl_find -- Find the index of a spectrum in the list, given its id
spl_get_spec -- Get a pointer to a spectrum from the list, given its id


ROUTINE DETAILS
Below, each routine is described in detail.

pointer spl_alloc (size)
Create a spectrum list

ARGUMENTS
size (input: integer)
Initial number of spectra the list will hold. If
unknown, specify INDEFI.

RETURNS
A pointer to the spectrum list structure. This pointer
should be used in all macros and routines which use the
spectrum list.

spl_free (o, free)
Destroy a spectrum list.

ARGUMENTS
o (input/output: pointer)
The pointer to the spectrum list structure to be
freed. On return, the value will be NULL.

free (input: bool)
If TRUE, each spectrum contained in the list will be
destroyed using the 'spunmap' routine. If false, the
spectra will not be destroyed. This is useful if there
are other structures which point to the same spectra.

spl_add (o, s, newid)
Add a spectrum to the list.

ARGUMENTS
o (input: pointer)
The spectrum list pointer.

s (input: pointer)
The pointer to the spectrum to add to the list.

newid (input: bool)
TRUE if a new id should be assigned to the component
being added.

spl_del (o, id, free)
Remove a spectrum from the list based on the spectrum id.

ARGUMENTS
o (input: pointer)
The spectrum list pointer.

id (input: int)
The id of the spectrum to remove from the list. If
INDEFI, the current spectrum will be removed. If there
is no current spectrum, an error will be generated.

free (input: bool)
If TRUE, the spectrum will be destroyed using the
'spunmap' routine. Else, only the list reference to
the spectrum will be destroyed, leaving the spectrum
itself intact. Useful if other memory structures refer
to the same spectrum.

pointer spl_copy (o, copy)
Create a new spectrum list that is a copy of the give one.

ARGUMENTS
o (input: pointer)
The spectrum list to be copied.

copy (input: bool)
If TRUE, a copy of each spectrum will be created using
the 'sp_copy' routine. The new list will then refer to
the copies of the spectra. Else, the new list will
contain the same pointers to the spectra. Hence, the
following scenario can occure: If a copy of a list is
made using 'copy = false', then the new list is
destroyed with 'free = true', the original list will
now have no spectra to point to.

RETURNS
A pointer to the new list is returned.

int spl_find (o, id)
Get the index into the array of a spectrum with the specified
id.

ARGUMENTS
o (input: pointer)
The spectrum list to search.

id (input: int)
The id of the spectrum to find. If INDEFI, the current
spectrum will be returned. If there is no current
spectrum, an error will be generated.

RETURNS
An index into the array pointing to the specified spectrum.

pointer spl_get_spec (o, id)
Get the spectrum from the list with the specified id.

ARGUMENTS
o (input: pointer)
The spectrum list to search.

id (input: int)
The id of the spectrum to retrieve. If INDEFI, the
current spectrum will be returned. If there is no
current spectrum, an error will be generated.

RETURNS
The pointer to the desired spectrum.

spl_merge (in, out)
Merge list IN with list OUT. This is based on id and version
number. The assumption is that list OUT is to contain the same
spectra as list IN (though attributes like units may be
different). "Sameness" is defined as similar id numbers and
similer SP_VER numbers. If an id exists in the OUT list that
doesn't exist in the IN list, it is delete. If spectra exist
in the IN list but not in the OUT list, they are added. For
spectra with the same SP_ID, if their SP_VER values are
different, the spectrum is deleted from the OUT list and
replaced with the one from the IN list.

ARGUMENTS
in [input: pointer to a spectrum list]
The spectrum list to merge into the spectrum list
specified in the out argument.

out [input/output: pointer to a spectrum list]
The spectrum list to be merged with. If NULL, a list
is created. The list can also be empty on input.



SPL_WAVE_DOMAIN (3Jan96) source SPL_WAVE_DOMAIN (3Jan96)



NAME
spl_wave_domain -- Find the wavelength domain of the list of spectra
sp_domain -- Find dispersion extrema of single spectrum in
specified units


DESCRIPTION
These routines determine the dispersion extrema of a spectrum or a
list of spectra. For spl_wave_domain, dispersion units of ANGSTROM
are assumed; for sp_domain the desired dispersion units are
specified on input, and the returned values are converted if
necessary. NOTE: "low" in this context means the minimum value of
the dispersion array in whatever units, and "high" means the
maximum value; whether these refer to the "blue" or the "red" end
of the spectrum depends upon the units of the spectrum or (for
sp_domain) the units specified in the calling argument.


ROUTINE DETAILS

procedure spl_wave_domain (s, low, high)


ARGUMENTS
s (input: pointer)
Spectrum List

low (output: double)
Lowest (shortest) wavelength

high (output: double)
Highest (longest) wavelength


procedure sp_domain (sp, dunit, low, high)


ARGUMENTS
sp (input: pointer)
Spectrum object

dunit (output: int)
index of desired units for low, high

low (output: double)
Lowest tabulated dispersion element

high (output: double)
Highest tabulated dispersion element



SPMAP (27Dec95) source SPMAP (27Dec95)



NAME
spmap -- allocate structure & arrays and initialize input spectrum


DESCRIPTION
Given a pointer to the parameter object:
1) Scan for a flux filename
2a) Try to open it as a table
or
2b) Try to open it as an image
3) Read the input spectrum and fill the spectrum object
Currently, we support reading:
SDAS binary tables containing flux and wavelength columns
ASCII tables containing flux and wavelength columns
STF multigroup format images, with individual images containing
flux and wavelength information
OIF images in the multispec format
FITS format OIF images in the multispec format
The routine reads all spectra that are avaible and returns a
pointer to a list of spectrum objects.


ROUTINE DETAILS

pointer procedure spmap (pp)


ARGUMENTS
pp (input: pointer)
Input pointer to parameter object

RETURNS
psl (output: pointer)
Output pointer to spectrum list



SP_MS_OPEN (27Dec95) source SP_MS_OPEN (27Dec95)



NAME
sp_ms_open -- Read a multispec format image into a spectrum object


DESCRIPTION
This routine takes an OIF image with WCS information and fills a
spectrum object. If the MWC is MULTISPEC, the routine will read
**all** of the spectra.


ROUTINE DETAILS

pointer procedure sp_ms_open ()


ARGUMENTS
ip (input: pointer)
input pointer to CIF object

im (input: pointer)
CIF primary file

RETURNS
A pointer to a spectrum list containing read spectra



SPNEAX (27Dec95) source SPNEAX (27Dec95)



NAME
sp_near_pix -- convert from wavelength to pixels for a spectrum


DESCRIPTION
This routine translates from wavelength back to pixel number in the
input spectrum.


ROUTINE DETAILS

int procedure sp_near_pix (s, w)


ARGUMENTS
s (input: pointer)
Spectrum object

w (input: double)
wavelength

RETURNS
An integer pixel value



SPPREP (27Dec95) source SPPREP (27Dec95)



NAME
sp_prep --


DESCRIPTION



ROUTINE DETAILS

procedure sp_prep (sp)
ARGUMENTS
sp (input: pointer)
Input spectrum pointer



SPTABN (27Dec95) SPTABN (27Dec95)



NAME
sp_table_open -- Read an ascii table or SDAS binary table and fill
a spectrum object


DESCRIPTION
This routine takes an table with columns containing flux,
wavelength, data quality and error values and fills a spectrum
object. The routine returns a pointer to a list of spectra (even
though only one spectrum is read in).


ROUTINE DETAILS

pointer procedure sp_table_open (filename, col_name, col_units)


ARGUMENTS
filename (input: char)
input filename

col_name (input: char)
array of default column names

col_units (input: char)
array of default column units

RETURNS
A pointer to a list of spectra



SPUNMAP (27Dec95) source SPUNMAP (27Dec95)



NAME
spunmap -- de-allocate dynamically allocated spectral arrays


DESCRIPTION
This routine destroys unwanted spectrum objects.


ROUTINE DETAILS

procedure spunmap (sp)


ARGUMENTS
sp (input: pointer)
Spectrum data structure



SPWRIE (1Feb96) source/specio SPWRIE (1Feb96)



NAME
spec_write -- write aspect spectrum object to disk file in STSDAS
table format


DESCRIPTION
This routine takes a spectrum (defined by a spectrum object in
memory) and writes it to disk as a SDAS binary table.


ROUTINE DETAILS

procedure sp_write (sp, outtab)


ARGUMENTS
sp (input: pointer)
Spectrum descriptor

outtab[SZ_FNAME] (input: char)
Output spectrum table name



TBINIT (27Dec95) source TBINIT (27Dec95)



NAME
tb_init -- initialize output STSDAS table columns


DESCRIPTION
This routine creates a table of "native" ASpec format for the
writing of spectra to disk.


ROUTINE DETAILS

procedure tb_init (tp, otcp)


ARGUMENTS
tp (input: pointer)
Table descriptor

otcp[NUM_COLS] (input: pointer)
Column descriptors



UPDATE_FILE_NAMES (27Dec95) source UPDATE_FILE_NAMES (27Dec95)



NAME
update_file_names -- Update filenames for use by the GUI


DESCRIPTION
This routine updates the parameter object with the names of files
read (while reading in a spectrum) to allow this information to be
sent back to the GUI.


ROUTINE DETAILS

procedure update_file_names (ip, pp)


ARGUMENTS
ip (input: pointer)
CIF pointer

pp (input: pointer)
Parameter object pointer



CL_GET_UNITS (Jan96) source/tasks CL_GET_UNITS (Jan96)



NAME
clg_unit_typ - Get enumerated type of units from CL string.
clg_units - Get units from CL parameters.


DESCRIPTION
clg_unit_typ --Interpret string to determine enumerated type of
units. Input string must be in uppercase. Unrecognized types
reported as zero.


ROUTINE DETAILS

procedure clg_unit_typ (unit_string)
ARGUMENTS
unit_string[ARB] (input: char)
string containing units and description

unit_type (output: int)
enumerated unti type; zero if unknown

procedure clg_units (dunit, funit)
ARGUMENTS
dunit, funit (output: int)
dispersion/flux unit enumerated type



T_AUTOSPEC (Jan96) source/tasks T_AUTOSPEC (Jan96)



NAME
t_autospec - Run the ASpect fitting in batch mode
get_cl_from_db - Build input component list from disk file & append to input CL
fit_report - Report fit progress to STDOUT


DESCRIPTION
Task to fit one or more model components to one or more input
spectra. Fetches task parameters, initializes data structures, etc.


ROUTINE DETAILS

procedure t_autospec ()
ARGUMENTS


procedure get_cl_from_db (cl, db_table)
ARGUMENTS
cl (input/output: pointer)
component list

db_table[ARB] (input: char)
name of input DB table of components


procedure fit_report (iter, chisq, force_flush)
ARGUMENTS
iter (input: int)
the current iteration count

chisq (input: real)
the current chisq value

force_flush (input: bool)
force flushing the output?



T_CVTSPEC (18Jan96) source/tasks T_CVTSPEC (18Jan96)



NAME
t_cvtspect -- Map FOS/HRS calibrated spectral images to TABLES
structure.


DESCRIPTION


ROUTINE DETAILS

procedure t_cvtspec ()
ARGUMENTS



T_EXPRESS (Jan96) source/tasks T_EXPRESS (Jan96)



NAME
t_express - Evaluate constraint expressions from a specified DB table
eval_expr - Evaluates component constraint expressions
exp_report - Prints component constraint expression evaluations


DESCRIPTION
Task to valuate constraint expressions from a specified DB table.
Takes input from disk tables, fetches task parameters, initializes
required data structures, etc.


ROUTINE DETAILS

procedure t_express ()
ARGUMENTS


procedure eval_expr (cl, funit, dunit, verbose)
ARGUMENTS
cl (input: pointer)
component list

dunit (input: int)
dispersion unit enumerated type

funit (input: int)
flux unit enumerated type

verbose (input: bool)
report in verbose mode?


procedure exp_report (cmp_index, cmp_name, par_name, val, expression, exp_status
ARGUMENTS
cmp_index (input: int)
component running index

cmp_name[ARB] (input: char)
component name

par_name[ARB] (input: char)
parameter name

val (input: double)
parameter value

expression[ARB] (input: char)
constraint expression

exp_status (input: short)
error status of expression evaluation



T_POPULATE (Jan96) source/tasks T_POPULATE (Jan96)



NAME
t_populate - Populate an ASpect DataBase entry with values from a p-set
fetch_pset - Open the correct p-set based on the component type
get_par_from_cl - Populate parameter values from a p-set


DESCRIPTION


ROUTINE DETAILS

procedure t_populate ()
ARGUMENTS


pointer procedure fetch_pset (comp_type)
ARGUMENTS
comp_type (input: int)
index of component type

pp (output: pointer)
pset descriptor for component type

RETURNS
pset name


procedure get_par_from_cl (pp, np, val, hi, lo, step, tol, sig, cnstr,
exprn)

ARGUMENTS
pp (input: pointer)
p-set descriptor

np (input: int)
number of parameters

val[ARB] (input: double)
parameter values

hi[ARB] (input: double)
upper bounds on parameter values

lo[ARB] (input: double)
lower bounds on parameter values

step[ARB] (input: real)
initial step size of parameter

tol[ARB] (input: real)
fit tolerance of parameter

sig[ARB] (input: real)
one-sigma errors on parameter

const[ARB] (input: short)
parameter constraint flag

exprn[SZ_EXPRESN+1, ARB] (input: char)
constraint flag



ABEDGE (Feb96) source/fitting/functions ABEDGE (Feb96)



NAME
abedge -- Absorption edge: optical depth ~(lambda/lambda_0)**3


DESCRIPTION
The routines in this library are for evaluating various physical or
empirical functions that describe components found in spectra of
astrophysical objects. The expectations for the units depend upon
the particular function (see list below): generally the physical
functions require specific units, while some others should be
defined by the user, and still others are independent of units, in
the sense that they can be defined in any unit system and should be
converted whenever the flux units for the evaluation changes.

Function Units Expected
-------- ----------------------------------------------------
abedge Expects dispersion in wavelengths (ANG would do)

The wavelengths are all assumed to be of type double precision
floating, and the fluxes and other arguments are single precision
floating.


ROUTINE DETAILS

procedure abedge (x, val, npts, x_refer, opdepth)
ARGUMENTS
x[ARB] (input: double)
evaluation wavelengths

val[ARB] (output: real)
array for evaluated function

npts (input: int)
number of points to evaluate

x_refer (input: double)
reference wavelength

opdepth (input: real)
optical depth normalization



ABGAUSS (Feb96) source/fitting/functions ABGAUSS (Feb96)



NAME
abgauss -- Gaussian absorption profile at wavelength X


DESCRIPTION
The routines in this library are for evaluating various physical or
empirical functions that describe components found in spectra of
astrophysical objects. The expectations for the units depend upon
the particular function (see list below): generally the physical
functions require specific units, while some others should be
defined by the user, and still others are independent of units, in
the sense that they can be defined in any unit system and should be
converted whenever the flux units for the evaluation changes.

Function Units Expected
-------- ----------------------------------------------------
abgauss Independent

The wavelengths are all assumed to be of type double precision
floating, and the fluxes and other arguments are single precision
floating.


ROUTINE DETAILS

procedure abgauss (x, val, npts, x_refer, ew, fwhm)
ARGUMENTS
x[ARB] (input: double)
evaluation wavelengths

val[ARB] (output: real)
array for evaluated function

npts (input: int)
number of points to evaluate

x_refer (input: double)
reference wavelength

ew (input: real)
equivalent width

fwhm (input: real)
full-width @half max



BPWRLAW (Mar96) source/fitting/functions BPWRLAW (Mar96)



NAME
bpwrlaw -- Broken (two-part) power-law continuum


DESCRIPTION
The routines in this library are for evaluating various physical or
empirical functions that describe components found in spectra of
astrophysical objects. The expectations for the units depend upon
the particular function (see list below): generally the physical
functions require specific units, while some others should be
defined by the user, and still others are independent of units, in
the sense that they can be defined in any unit system and should be
converted whenever the flux units for the evaluation changes.

Function Units Expected
-------- ----------------------------------------------------
bpwrlaw User defined

The wavelengths are all assumed to be of type double precision
floating, and the fluxes and other arguments are single precision
floating.


ROUTINE DETAILS

procedure bpwrlaw (x, val, npts, x_refer, flux, alpha1, alpha2)
ARGUMENTS
x[ARB] (input: double)
evaluation wavelengths

val[ARB] (output: real)
array for evaluated function

.ls x_refer (input: double)
reference wavelength

npts (input: int)
number of points to evaluate

flux (input: real)
flux @ref wavelength

alpha1 (input: real)
powerlaw indexes above "x_refer"

alpha2 (input: real)
powerlaw indexes below "x_refer"



DAMPABS (Feb96) source/fitting/functions DAMPABS (Feb96)



NAME
dampabs -- Damped absorption profile


DESCRIPTION
The routines in this library are for evaluating various physical or
empirical functions that describe components found in spectra of
astrophysical objects. The expectations for the units depend upon
the particular function (see list below): generally the physical
functions require specific units, while some others should be
defined by the user, and still others are independent of units, in
the sense that they can be defined in any unit system and should be
converted whenever the flux units for the evaluation changes.

Function Units Expected
-------- ----------------------------------------------------
dampabs Expects dispersion in Angstroms (but converts to Hz)

The wavelengths are all assumed to be of type double precision
floating, and the fluxes and other arguments are single precision
floating.


ROUTINE DETAILS

procedure dampabs (x, val, npts, x_refer, nhf, gamma)
ARGUMENTS
,ls x[ARB] (input: double) evaluation wavelength array

val[ARB] (output: real)
array for evaluated function

npts (input: int)
number of points to evaluate

nhf (input: double)
product of column density & oscillator strength

gamma (input: real)
Transition lifetime (sec)



EXPONEN (Feb96) source/fitting/functions EXPONEN (Feb96)



NAME
exponen -- Exponential profile


DESCRIPTION
The routines in this library are for evaluating various physical or
empirical functions that describe components found in spectra of
astrophysical objects. The expectations for the units depend upon
the particular function (see list below): generally the physical
functions require specific units, while some others should be
defined by the user, and still others are independent of units, in
the sense that they can be defined in any unit system and should be
converted whenever the flux units for the evaluation changes.

Function Units Expected
-------- ----------------------------------------------------
exponen Independent

The wavelengths are all assumed to be of type double precision
floating, and the fluxes and other arguments are single precision
floating.


ROUTINE DETAILS

procedure exponen (x, val, npts, x_refer, ew, fwhm)
ARGUMENTS
x[ARB] (input: double)
evaluation wavelengths

val[ARB] (output: real)
array for evaluated function

npts (input: int)
number of points to evaluate

x_refer (input: double)
reference wavelength

ew (input: real)
equivalent width

fwhm (input: real)
full width-half max



EXTINCT (Feb96) source/fitting/functions EXTINCT (Feb96)



NAME
extinct -- Interstellar extinction correction @wavelength X


DESCRIPTION
The routines in this library are for evaluating various physical or
empirical functions that describe components found in spectra of
astrophysical objects. The expectations for the units depend upon
the particular function (see list below): generally the physical
functions require specific units, while some others should be
defined by the user, and still others are independent of units, in
the sense that they can be defined in any unit system and should be
converted whenever the flux units for the evaluation changes.

Function Units Expected
-------- ----------------------------------------------------
extinct Expects dispersion in Angstroms

The wavelengths are all assumed to be of type double precision
floating, and the fluxes and other arguments are single precision
floating.


ROUTINE DETAILS

procedure extinct (x, val, npts, ebmv, redfunc, rv)
ARGUMENTS
x[ARB] (input: double)
evaluation wavelength array

val[ARB] (output: real)
array for evaluated function

npts (input: int)
number of points to evaluate

ebmv (input: real)
E(B-V): color excess

redfunc (input: int)
choice of redding law

rv (input: real)
ratio of total to selective extinction (@V)



GAUSS (Feb96) source/fitting/functions GAUSS (Feb96)



NAME
gauss -- Gaussian emission profile @wavelength X


DESCRIPTION
The routines in this library are for evaluating various physical or
empirical functions that describe components found in spectra of
astrophysical objects. The expectations for the units depend upon
the particular function (see list below): generally the physical
functions require specific units, while some others should be
defined by the user, and still others are independent of units, in
the sense that they can be defined in any unit system and should be
converted whenever the flux units for the evaluation changes.

Function Units Expected
-------- ----------------------------------------------------
gauss Independent

The wavelengths are all assumed to be of type double precision
floating, and the fluxes and other arguments are single precision
floating.


ROUTINE DETAILS

procedure gauss (x, val, npts, x_refer, flux, fwhm, skew)
ARGUMENTS
x[ARB] (input: double)
evaluation wavelengths

val[ARB] (output: real)
array for evaluated function

npts (input: int)
number of points to evaluate

x_refer (input: double)
reference wavelength

flux (input: real)
integrated flux

fwhmw (input: real)
full width in km/s

skew (input: real)
skew



LOGABS (Feb96) source/fitting/functions LOGABS (Feb96)



NAME
logabs -- Absorption profile of the form: fmax * (x/x0)**(-alpha)


DESCRIPTION
The routines in this library are for evaluating various physical or
empirical functions that describe components found in spectra of
astrophysical objects. The expectations for the units depend upon
the particular function (see list below): generally the physical
functions require specific units, while some others should be
defined by the user, and still others are independent of units, in
the sense that they can be defined in any unit system and should be
converted whenever the flux units for the evaluation changes.

Function Units Expected
-------- ----------------------------------------------------
logabs Expects dispersion in wavelengths (ANG would do)

The wavelengths are all assumed to be of type double precision
floating, and the fluxes and other arguments are single precision
floating.


ROUTINE DETAILS

procedure logabs (x, val, npts, x_refer, flux, sigma)
ARGUMENTS
x[ARB] (input: double)
evaluation wavelengths

val[ARB] (output: real)
array for evaluated function

npts (input: int)
number of points to evaluate

x_refer (input: double)
reference wavelength

flux (input: real)
flux normalization

sigma (input: real)
feature width



LOGARITH (Feb96) source/fitting/functions LOGARITH (Feb96)



NAME
logarith -- Power-law line profile: F = I_o (x/centroid)**alpha


DESCRIPTION
The routines in this library are for evaluating various physical or
empirical functions that describe components found in spectra of
astrophysical objects. The expectations for the units depend upon
the particular function (see list below): generally the physical
functions require specific units, while some others should be
defined by the user, and still others are independent of units, in
the sense that they can be defined in any unit system and should be
converted whenever the flux units for the evaluation changes.

Function Units Expected
-------- ----------------------------------------------------
logarith Expects dispersion in wavelengths (ANG would do)

The wavelengths are all assumed to be of type double precision
floating, and the fluxes and other arguments are single precision
floating.


ROUTINE DETAILS

procedure logarith (x, val, npts, x_refer, flux, sigma, skew)
ARGUMENTS
x[ARB] (input: double)
evaluation wavelengths

val[ARB] (output: real)
array for evaluated function

npts (input: int)
number of points to evaluate

x_refer (input: double)
reference wavelength

flux (input: real)
flux normalization

sigma (input: real)
feature width

skew (input: real)
skew



LORENTZ (Feb96) source/fitting/functions LORENTZ (Feb96)



NAME
lorentz -- Modified Lorentzian profile


DESCRIPTION
The routines in this library are for evaluating various physical or
empirical functions that describe components found in spectra of
astrophysical objects. The expectations for the units depend upon
the particular function (see list below): generally the physical
functions require specific units, while some others should be
defined by the user, and still others are independent of units, in
the sense that they can be defined in any unit system and should be
converted whenever the flux units for the evaluation changes.

Function Units Expected
-------- ----------------------------------------------------
lorentz Independent

The wavelengths are all assumed to be of type double precision
floating, and the fluxes and other arguments are single precision
floating.


ROUTINE DETAILS

procedure lorentz (x, val, npts, x_refer, norm, fwhm, alpha)
ARGUMENTS
x[ARB] (input: double)
evaluation wavelengths

val[ARB] (output: real)
array for evaluated function

npts (input: int)
number of points to evaluate

x_refer (input: double)
reference wavelength

norm (input: real)
flux normalization

fwhm (input: real)
feature width

alpha (input: real)
power-law index.



PLANCK (Feb96) source/fitting/functions PLANCK (Feb96)



NAME
planck -- Black-body continuum


DESCRIPTION
The routines in this library are for evaluating various physical or
empirical functions that describe components found in spectra of
astrophysical objects. The expectations for the units depend upon
the particular function (see list below): generally the physical
functions require specific units, while some others should be
defined by the user, and still others are independent of units, in
the sense that they can be defined in any unit system and should be
converted whenever the flux units for the evaluation changes.

Function Units Expected
-------- ----------------------------------------------------
planck Expects wavelength in ANG,
flux defined per unit wavelength

The wavelengths are all assumed to be of type double precision
floating, and the fluxes and other arguments are single precision
floating.


ROUTINE DETAILS

procedure planck (x, val, npts, x_refer, flux, temperature)
ARGUMENTS
x[ARB] (input: double)
evaluation wavelengths

val[ARB] (output: real)
array for evaluated function

npts (input: int)
number of points to evaluate

x_refer (input: double)
reference wavelength

flux (input: real)
flux @reference wavelength

temperature (input: real)
temperature (Kelvin)



POLY (Feb96) source/fitting/functions POLY (Feb96)



NAME
poly -- Polynomial continuum


DESCRIPTION
The routines in this library are for evaluating various physical or
empirical functions that describe components found in spectra of
astrophysical objects. The expectations for the units depend upon
the particular function (see list below): generally the physical
functions require specific units, while some others should be
defined by the user, and still others are independent of units, in
the sense that they can be defined in any unit system and should be
converted whenever the flux units for the evaluation changes.

Function Units Expected
-------- ----------------------------------------------------
poly User defined

The wavelengths are all assumed to be of type double precision
floating, and the fluxes and other arguments are single precision
floating.


ROUTINE DETAILS

procedure poly (x, val, npts, x_refer, coef, order)
ARGUMENTS
x[ARB] (input: double)
evaluation wavelengths

val[ARB] (output: real)
array for evaluated function

npts (input: int)
number of points to evaluate

x_refer (input: double)
reference wavelength

coef[ARB] (input: real)
polynomial coefficients 0 through ORDER

order (input: int)
order of polynomial (i.e, highest power of X)



PWRLAW (Feb96) source/fitting/functions PWRLAW (Feb96)



NAME
pwrlaw -- Power-law continuum


DESCRIPTION
The routines in this library are for evaluating various physical or
empirical functions that describe components found in spectra of
astrophysical objects. The expectations for the units depend upon
the particular function (see list below): generally the physical
functions require specific units, while some others should be
defined by the user, and still others are independent of units, in
the sense that they can be defined in any unit system and should be
converted whenever the flux units for the evaluation changes.

Function Units Expected
-------- ----------------------------------------------------
pwrlaw User defined

The wavelengths are all assumed to be of type double precision
floating, and the fluxes and other arguments are single precision
floating.


ROUTINE DETAILS

procedure pwrlaw (x, val, npts, x_refer, flux, alpha)
ARGUMENTS
x[ARB] (input: double)
evaluation wavelengths

val[ARB] (output: real)
array for evaluated function

npts (input: int)
number of points to evaluate

x_refer (input: double)
reference wavelength

flux (input: real)
flux @ref wave

alpha (input: real)
spectral index



RECOMB (Feb96) source/fitting/functions RECOMB (Feb96)



NAME
recomb -- Hydrogen Recombination continuum


DESCRIPTION
The routines in this library are for evaluating various physical or
empirical functions that describe components found in spectra of
astrophysical objects. The expectations for the units depend upon
the particular function (see list below): generally the physical
functions require specific units, while some others should be
defined by the user, and still others are independent of units, in
the sense that they can be defined in any unit system and should be
converted whenever the flux units for the evaluation changes.

Function Units Expected
-------- ----------------------------------------------------
recomb Expects wavelength in ANG,

The wavelengths are all assumed to be of type double precision
floating, and the fluxes and other arguments are single precision
floating.


ROUTINE DETAILS

procedure recomb (x, val, npts, x_refer, flux, temperature, fwhm)
ARGUMENTS
x[ARB] (input: double)
evaluation wavelengths

val[ARB] (output: real)
array for evaluated function

npts (input: int)
number of points to evaluate

x_refer (input: double)
reference wavelength

flux (input: real)
flux at reference wavelength

temperature (input: real)
temperature (Kelvin)

fwhm (input: real)



REDLAW (Jan96) nebular/lib REDLAW (Jan96)



NAME
gal_redlaw -- Calculates the ISEF for the Galaxy by Savage & Mathis
(1979) lmc_redlaw -- Calculates the ISEF for the LMC by Howarth
(1983) smc_redlaw -- Calculates the ISEF for the SMC by Prevot, et
al. (1984) ccm_redlaw -- Calculates the ISEF for the Galaxy by
Cardelli, Clayton & Mathis (1989)


DESCRIPTION
The following procedures return an evaluation of the interstellar
extinction function (ISEF) as published by various authors. Often
the ISEF is expressed as A(lambda)/E(B-V), and is normalized such
that the function is zero at the V passband, and 1.00 at B. These
functions are expressed in the literature as tabulated functions,
polynomial expressions, or both. The parameter is wavelength in
inverse microns, and the result is color excess in magnitudes. The
functions are defined or extended to about lambda=1000 Ang; the
form of the ISEF is not known shortward of that.

The routines here return A(lambda)/A(V) by renormaling the
extinction function to 0.0 at zero inverse wavelength, and 1.0 at
V. The ratio of total to selective absorption in the V passband,
R_v, is thought to average 3.10 for most lines of sight in the
Galaxy and in the LMC; it may be somewhat lower in the SMC. It is
the responsibility of the calling function to re-normalize (i.e.,
multiply) the returned value by the choice of R_v. Note that R_v
is not a parameter for any of the ISEF's here except for CCM.


ROUTINE DETAILS

procedure gal_redlaw (wave, extl, npts)
ARGUMENTS
wave[ARB] (input: double)
wavelength of interest

extl[ARB] (output: real)
extinction evaluation array

npts (input: int)
size of evaluation/wave arrays


procedure ccm_redlaw (wave, extl, npts, rv)
ARGUMENTS
wave[ARB] (input: double)
wavelength of interest

extl[ARB] (output: real)
extinction evaluation array

npts (input: int)
size of evaluation/wave arrays

rv (input: real)
ratio of total to selective extinction


procedure lmc_redlaw (wave, extl, npts)
ARGUMENTS
wave[ARB] (input: double)
wavelength of interest

extl[ARB] (output: real)
extinction evaluation array

npts (input: int)
size of evaluation/wave arrays


procedure smc_redlaw (wave, extl, npts)
ARGUMENTS
wave[ARB] (input: double)
wavelength of interest

extl[ARB] (output: real)
extinction evaluation array

npts (input: int)
size of evaluation/wave arrays



USR_ABS (Feb96) source/fitting/functions USR_ABS (Feb96)



NAME
usr_abs -- Evaluate a user-defined (tabulated) absorption feature.


DESCRIPTION
The routines in this library are for evaluating various physical or
empirical functions that describe components found in spectra of
astrophysical objects. The expectations for the units depend upon
the particular function (see list below): generally the physical
functions require specific units, while some others should be
defined by the user, and still others are independent of units, in
the sense that they can be defined in any unit system and should be
converted whenever the flux units for the evaluation changes.

Function Units Expected
-------- ----------------------------------------------------
usr_abs User defined

The wavelengths are all assumed to be of type double precision
floating, and the fluxes and other arguments are single precision
floating.


ROUTINE DETAILS

procedure usr_abs (x, val, npts, usr_x, usr_func, sz_tab, norm, shift, z)
ARGUMENTS
x[ARB] (input: double)
evaluation wavelength array

val[ARB (output: real)
array for evaluated function

npts (input: int)
number of points to evaluate

usr_x[ARB] (input: double)

usr_func[ARB] (input: real)
tabulated function array

sz_tab (input: int)
size of tabulated function array

norm (input: real)
flux normalization

shift (input: double)
wavelength shift

z (input: double)
wavelength stretch (redshift)



USR_EMISS (Feb96) source/fitting/functions USR_EMISS (Feb96)



NAME
usr_emiss -- Evaluate a user-defined (tabulated) emission feature.


DESCRIPTION
The routines in this library are for evaluating various physical or
empirical functions that describe components found in spectra of
astrophysical objects. The expectations for the units depend upon
the particular function (see list below): generally the physical
functions require specific units, while some others should be
defined by the user, and still others are independent of units, in
the sense that they can be defined in any unit system and should be
converted whenever the flux units for the evaluation changes.

Function Units Expected
-------- ----------------------------------------------------
usr_emiss User defined

The wavelengths are all assumed to be of type double precision
floating, and the fluxes and other arguments are single precision
floating.


ROUTINE DETAILS

procedure usr_emiss (x, val, npts, usr_x, usr_func, sz_tab, norm, shift, z)
ARGUMENTS
x[ARB] (input: double)
evaluation wavelength array

val[ARB (output: real)
array for evaluated function

npts (input: int)
number of points to evaluate

usr_x[ARB] (input: double)

usr_func[ARB] (input: real)
tabulated function array

sz_tab (input: int)
size of tabulated function array

norm (input: real)
flux normalization

shift (input: double)
wavelength shift

z (input: double)
wavelength stretch (redshift)



VOIGT (Feb96) source/fitting/functions VOIGT (Feb96)



NAME
voigt -- Evaluate a Voigt line profile (from Lang 1980)


DESCRIPTION
The routines in this library are for evaluating various physical or
empirical functions that describe components found in spectra of
astrophysical objects. The expectations for the units depend upon
the particular function (see list below): generally the physical
functions require specific units, while some others should be
defined by the user, and still others are independent of units, in
the sense that they can be defined in any unit system and should be
converted whenever the flux units for the evaluation changes.

Function Units Expected
-------- ----------------------------------------------------
voigt Expects wavelength in ANG (but converts to Hz)

The wavelengths are all assumed to be of type double precision
floating, and the fluxes and other arguments are single precision
floating.


ROUTINE DETAILS

procedure voigt (x, val, npts, x_refer, flux, width_d, width_l)
ARGUMENTS
x[ARB] (input: double)
evaluation wavelength array

val[ARB] (output: real)
array for evaluated function

npts (input: int)
number of points to evaluate

x_refer (input: double)
reference wavelength

flux (input: real)
line flux

width_d (input: real)
Doppler widths in km/s

width_l (input: real)
Lorentz widths in km/s



NAME
voigt_profile - Compute the real portion of the Voigt profile.
#@(#) funits.x 1.1
include
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, nscan(), strdic()
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)
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
#@(#) getimage.x 1.1
include
include
include "smw.h"

# GETIMAGE -- Read new image pixels.

procedure getimage (image, nline, nband, nap, wave_scl, w0, wpc, units,
im, mw, sh)

char image[ARB]
int nline, nband, nap
bool wave_scl
double w0, wpc
char units[ARB]
pointer sp, imsect, im, mw, sh

int i, n, sec[3,3], clgeti()
pointer immap(), smw_openim()
errchk immap, shdr_open, shdr_system, un_changer

begin
call smark (sp)
call salloc (imsect, SZ_FNAME, TY_CHAR)

# Map the image if necessary. Don't allow image sections but
# determine requested spectrum from any explicit specification.

if (im == NULL) {
call imgsection (image, Memc[imsect], SZ_FNAME)
call imgimage (image, image, SZ_FNAME)
im = immap (image, READ_ONLY, 0)
mw = smw_openim (im)
n = IM_NDIM(im)
if (Memc[imsect] != EOS && n > 1) {
call amovki (1, sec[1,1], n)
call amovki (IM_LEN(im,1), sec[1,2], n)
call amovki (1, sec[1,3], n)
call id_section (Memc[imsect], sec[1,1], sec[1,2], sec[1,3], n)
switch (SMW_FORMAT(mw)) {
case SMW_ND:
i = 0
if (n == 2) {
if (abs (sec[1,2]-sec[1,1]) == 0) {
nline = sec[1,1]
i = 2
} else if (abs (sec[2,2]-sec[2,1]) == 0) {
nline = sec[2,1]
i = 1
}
} else {
if (abs (sec[1,2]-sec[1,1]) == 0) {
nline = sec[1,1]
if (abs (sec[2,2]-sec[2,1]) == 0) {
nband = sec[2,1]
if (abs (sec[3,2]-sec[3,1]) > 0)
i = 3
} else if (abs (sec[3,2]-sec[3,1]) == 0) {
nband = sec[3,1]
i = 2
}
} else if (abs (sec[2,2]-sec[2,1]) == 0) {
nline = sec[2,1]
if (abs (sec[3,2]-sec[3,1]) == 0) {
nband = sec[3,1]
i = 1
}
}
}
if (i > 0) {
call smw_daxis (mw, im, i, INDEFI, INDEFI)
call smw_saxis (mw, NULL, im)
}
default:
if (abs (sec[2,2]-sec[2,1]) == 0)
nline = sec[2,1]
if (n > 2 && abs (sec[3,2]-sec[3,1]) == 0)
nband = sec[3,1]
}
}
}

# Get header info.
switch (SMW_FORMAT(mw)) {
case SMW_ND:
nap = INDEFI
n = SMW_LLEN(mw,2)
if (n > 1) {
if (nline == 0)
nline = max (1, min (n, clgeti ("line")))
} else
nline = 0
n = SMW_LLEN(mw,3)
if (n > 1) {
if (nband == 0)
nband = max (1, min (n, clgeti ("band")))
} else
nband = 0
default:
n = SMW_NSPEC(mw)
if (n > 1) {
if (nline == 0) {
nline = clgeti ("line")
nap = nline
}
} else {
nline = 0
nap = INDEFI
}
n = SMW_NBANDS(mw)
if (n > 1) {
if (nband == 0)
nband = max (1, min (n, clgeti ("band")))
} else
nband = 0
}

call shdr_open (im, mw, nline, nband, nap, SHDATA, sh)
nap = AP(sh)

if (DC(sh) == DCNO && !IS_INDEFD(w0))
call usercoord (sh, 'l', 1D0, w0, 2D0, w0+wpc)

# Cancel wavelength coordinates if not desired or set units.
if (!wave_scl)
call shdr_system (sh, "physical")
else {
iferr (call shdr_units (sh, units))
;
}

call sfree (sp)
end
#@(#) idmap.x 1.1
include
include
include
include "smw.h"
include "identify.h"

# Sepcial section words.
define SPECIAL "|first|middle|x|y|z|last|column|line|band|"
define FIRST 1
define MIDDLE 2
define X 3
define Y 4
define Z 5
define LAST 6
define COLUMN 7
define LINE 8
define BAND 9

# ID_MAP -- Map an image for IDENTIFY/REIDENTIFY
# The image must 1, 2, or 3 dimensional. An image section may be given with
# the image name or with the CL parameter "section". The CL parameter can
# have one of the following formats:
# 1. An IMIO image section
# 2. [line|column|x|y|z] [#|middle|last] [#|middle|last]
# 3. [#|middle|last] [#|middle|last] [line|column|x|y|z]
# where # is a line or column number. The strings may be abbreviated.
# The task returns and error if it cannot map the image or determine
# the 1D line or column desired.

procedure id_map (id)

pointer id # IDENTIFY data structure pointer

int i, j, k, l, a, b, c, x1[3], x2[3], xs[3]
pointer sp, wrd1, wrd2, wrd3, im

int imaccess(), strdic(), ctoi(), nscan()
pointer immap()
errchk immap, id_maphdr

begin
# Separate the image name and image section and map the full image.
call imgsection (Memc[ID_IMAGE(id)], Memc[ID_SECTION(id)], SZ_FNAME)
call imgimage (Memc[ID_IMAGE(id)], Memc[ID_IMAGE(id)], SZ_FNAME)
im = immap (Memc[ID_IMAGE(id)], READ_ONLY, 0)

# If no image section is found use the "section" parameter.
if (Memc[ID_SECTION(id)] == EOS && IM_NDIM(im) > 1) {
call clgstr ("section", Memc[ID_SECTION(id)], SZ_FNAME)
call xt_stripwhite (Memc[ID_SECTION(id)])

# If not an image section construct one.
if (Memc[ID_SECTION(id)] != '[') {
call smark (sp)
call salloc (wrd1, SZ_FNAME, TY_CHAR)
call salloc (wrd2, SZ_FNAME, TY_CHAR)
call salloc (wrd3, SZ_FNAME, TY_CHAR)

call sscan (Memc[ID_SECTION(id)])

# Parse axis and elements.
call gargwrd (Memc[wrd1], SZ_FNAME)
call gargwrd (Memc[wrd2], SZ_FNAME)
call gargwrd (Memc[wrd3], SZ_FNAME)
switch (nscan()) {
case 0:
a = X
b = MIDDLE
c = MIDDLE
case 1:
a = strdic (Memc[wrd1], Memc[wrd1], SZ_FNAME, SPECIAL)
b = MIDDLE
c = MIDDLE
case 2:
a = strdic (Memc[wrd1], Memc[wrd1], SZ_FNAME, SPECIAL)
if (a >= X)
b = strdic (Memc[wrd2], Memc[wrd2], SZ_FNAME, SPECIAL)
else {
b = a
a = strdic (Memc[wrd2], Memc[wrd2], SZ_FNAME, SPECIAL)
call strcpy (Memc[wrd1], Memc[wrd2], SZ_FNAME)
}
c = MIDDLE
call strcpy (Memc[wrd2], Memc[wrd3], SZ_FNAME)
case 3:
a = strdic (Memc[wrd1], Memc[wrd1], SZ_FNAME, SPECIAL)
if (a >= X) {
b = strdic (Memc[wrd2], Memc[wrd2], SZ_FNAME, SPECIAL)
c = strdic (Memc[wrd3], Memc[wrd3], SZ_FNAME, SPECIAL)
} else {
b = a
a = strdic (Memc[wrd2], Memc[wrd2], SZ_FNAME, SPECIAL)
if (a >= X) {
c = strdic (Memc[wrd3], Memc[wrd3],SZ_FNAME,SPECIAL)
call strcpy (Memc[wrd1], Memc[wrd2], SZ_FNAME)
} else {
c = b
b = a
a = strdic (Memc[wrd3], Memc[wrd3],SZ_FNAME,SPECIAL)
call strcpy (Memc[wrd2], Memc[wrd3], SZ_FNAME)
call strcpy (Memc[wrd1], Memc[wrd2], SZ_FNAME)
}
}
}

switch (a) {
case X, LINE:
i = 1
j = 2
k = 3
case Y, COLUMN:
i = 2
j = 1
k = 3
case Z, BAND:
i = 3
j = 1
k = 2
default:
call imunmap (im)
call error (1,
"Error in section specification or non-unique abbreviation")
}

switch (b) {
case FIRST:
ID_LINE(id,1) = 1
case MIDDLE:
ID_LINE(id,1) = (1 + IM_LEN(im,j)) / 2
case LAST:
ID_LINE(id,1) = IM_LEN(im,j)
default:
l = 1
if (ctoi (Memc[wrd2], l, ID_LINE(id,1)) == 0)
call error (1, "Error in section specification")
}

switch (c) {
case FIRST:
ID_LINE(id,2) = 1
case MIDDLE:
ID_LINE(id,2) = (1 + IM_LEN(im,k)) / 2
case LAST:
ID_LINE(id,2) = IM_LEN(im,k)
default:
l = 1
if (ctoi (Memc[wrd3], l, ID_LINE(id,2)) == 0)
call error (1, "Error in section specification")
}

# Format section.
switch (IM_NDIM(im)) {
case 2:
switch (i) {
case 1:
call sprintf (Memc[ID_SECTION(id)], SZ_FNAME, "[*,%d]")
case 2:
call sprintf (Memc[ID_SECTION(id)], SZ_FNAME, "[%d,*]")
default:
call error (1, "Error in section specification")
}
call pargi (ID_LINE(id,1))
case 3:
switch (i) {
case 1:
call sprintf (Memc[ID_SECTION(id)],SZ_FNAME,"[*,%d,%d]")
case 2:
call sprintf (Memc[ID_SECTION(id)],SZ_FNAME,"[%d,*,%d]")
case 3:
call sprintf (Memc[ID_SECTION(id)],SZ_FNAME,"[%d,%d,*]")
}
call pargi (ID_LINE(id,1))
call pargi (ID_LINE(id,2))
case 4:
call error (1, "Image dimension greater than 3 not allowed")
}
}
}

# Parse the image section.
x1[1] = 1; x2[1] = IM_LEN(im,1); xs[1] = 1
x1[2] = 1; x2[2] = IM_LEN(im,2); xs[2] = 1
x1[3] = 1; x2[3] = IM_LEN(im,3); xs[3] = 1
call id_section (Memc[ID_SECTION(id)], x1, x2, xs, 3)

# Set the axes. The axis to be identified is the longest one.
i = 1
if (IM_NDIM(im) > 1 && abs (x1[2]-x2[2]) >= abs (x1[i]-x2[i]))
i = 2
if (IM_NDIM(im) > 2 && abs (x1[3]-x2[3]) >= abs (x1[i]-x2[i]))
i = 3
if (IM_NDIM(im) > 3)
call error (1, "Image dimension greater than 3 not allowed")

switch (i) {
case 1:
j = 2
k = 3
case 2:
j = 1
k = 3
case 3:
j = 1
k = 2
}

ID_LINE(id,1) = (x1[j] + x2[j]) / 2
ID_LINE(id,2) = (x1[k] + x2[k]) / 2
ID_MAXLINE(id,1) = IM_LEN(im, j)
ID_MAXLINE(id,2) = IM_LEN(im, k)
ID_NSUM(id,1) = min (ID_MAXLINE(id,1), ID_NSUM(id,1))
ID_NSUM(id,2) = min (ID_MAXLINE(id,2), ID_NSUM(id,2))
call smw_daxis (NULL, NULL, i, ID_NSUM(id,1), ID_NSUM(id,2))

call id_maphdr (id, im)

# Open the image READ_WRITE if possible in order to add REFSPEC.
# This is not done earlier to avoid updating of the WCS.

call imunmap (im)
if (imaccess (Memc[ID_IMAGE(id)], READ_WRITE) == YES)
im = immap (Memc[ID_IMAGE(id)], READ_WRITE, 0)
else
im = immap (Memc[ID_IMAGE(id)], READ_ONLY, 0)
call id_noextn (Memc[ID_IMAGE(id)])
IM(ID_SH(id)) = im
end


# ID_MAPHDR -- Map image header.

procedure id_maphdr (id, im)

pointer id # ID pointer
pointer im # IMIO pointer

int i
pointer mw, sh, smw_openim(), smw_sctran()
errchk smw_openim(), shdr_open(), smw_sctran

begin
mw = smw_openim (im)
if (SMW_TRANS(mw) == YES) {
if (SMW_PAXIS(mw,1) == 1)
call smw_daxis (mw, im, 2, INDEFI, INDEFI)
else
call smw_daxis (mw, im, 1, INDEFI, INDEFI)
call smw_saxes (mw, NULL, im)
}
call shdr_open (im, mw, ID_LINE(id,1), ID_LINE(id,2),
INDEFI, SHHDR, ID_SH(id))
sh = ID_SH(id)

if (SMW_FORMAT(mw) == SMW_MS || SMW_FORMAT(mw) == SMW_ES) {
ID_MAXLINE(id,1) = IM_LEN(im,2)
ID_MAXLINE(id,2) = IM_LEN(im,3)
ID_NSUM(id,1) = 1
ID_NSUM(id,2) = 1
ID_LINE(id,1) = max (1, min (ID_MAXLINE(id,1), ID_LINE(id,1)))
ID_LINE(id,2) = 1
call mfree (ID_APS(id), TY_INT)
call malloc (ID_APS(id), ID_MAXLINE(id,1), TY_INT)
do i = 1, ID_MAXLINE(id,1) {
call shdr_open (im, mw, i, 1, INDEFI, SHHDR, sh)
Memi[ID_APS(id)+i-1] = AP(sh)
}
ID_AP(id,1) = Memi[ID_APS(id)+ID_LINE(id,1)-1]
ID_AP(id,2) = 1
} else {
call mfree (ID_APS(id), TY_INT)
ID_AP(id,1) = ID_LINE(id,1)
ID_AP(id,2) = ID_LINE(id,2)
}
ID_NPTS(id) = IM_LEN(im, SMW_LAXIS(mw,1))

call gt_sets (ID_GT(id), GTXLABEL, LABEL(sh))
call ic_pstr (ID_IC(id), "ylabel", LABEL(sh))
call gt_sets (ID_GT(id), GTXUNITS, UNITS(sh))
call ic_pstr (ID_IC(id), "yunits", UNITS(sh))

# Set logical / physical transformations
i = 2 ** (SMW_PAXIS(mw,1) - 1)
ID_LP(id) = smw_sctran (mw, "logical", "physical", i)
ID_PL(id) = smw_sctran (mw, "physical", "logical", i)
end


# ID_SECTION -- Parse an image section into its elements.
# 1. The default values must be set by the caller.
# 2. A null image section is OK.
# 3. The first nonwhitespace character must be '['.
# 4. The last interpreted character must be ']'.
#
# This procedure should be replaced with an IMIO procedure at some
# point.

procedure id_section (section, x1, x2, xs, ndim)

char section[ARB] # Image section
int x1[ndim] # Starting pixel
int x2[ndim] # Ending pixel
int xs[ndim] # Step
int ndim # Number of dimensions

int i, ip, a, b, c, temp, ctoi()
define error_ 99

begin
# Decode the section string.
ip = 1
while (IS_WHITE(section[ip]))
ip = ip + 1
if (section[ip] == '[')
ip = ip + 1
else if (section[ip] == EOS)
return
else
goto error_

do i = 1, ndim {
while (IS_WHITE(section[ip]))
ip = ip + 1
if (section[ip] == ']')
break

# Default values
a = x1[i]
b = x2[i]
c = xs[i]

# Get a:b:c. Allow notation such as "-*:c"
# (or even "-:c") where the step is obviously negative.

if (ctoi (section, ip, temp) > 0) { # a
a = temp
if (section[ip] == ':') {
ip = ip + 1
if (ctoi (section, ip, b) == 0) # a:b
goto error_
} else
b = a
} else if (section[ip] == '-') { # -*
temp = a
a = b
b = temp
ip = ip + 1
if (section[ip] == '*')
ip = ip + 1
} else if (section[ip] == '*') # *
ip = ip + 1
if (section[ip] == ':') { # ..:step
ip = ip + 1
if (ctoi (section, ip, c) == 0)
goto error_
else if (c == 0)
goto error_
}
if (a > b && c > 0)
c = -c

x1[i] = a
x2[i] = b
xs[i] = c

while (IS_WHITE(section[ip]))
ip = ip + 1
if (section[ip] == ',')
ip = ip + 1
}

if (section[ip] != ']')
goto error_

return
error_
call error (0, "Error in image section specification")
end
#@(#) idnoextn.x 1.1
# ID_NOEXTN -- Remove standard image extensions.

procedure id_noextn (image)

char image[ARB] # Image name

int strlen()

begin
call xt_imroot (image, image, strlen (image))
end
#@(#) mktitle.x 1.1
include
include "smw.h"
include "units.h"

# MKTITLE -- Make a spectrum title (IIDS style)

procedure mktitle (sh, gt)

pointer sh, gt

pointer sp, str

begin
# Do nothing if the GTOOLS pointer is undefined.
if (gt == NULL)
return

call smark (sp)
call salloc (str, SZ_LINE, TY_CHAR)

call sprintf (Memc[str], SZ_LINE,
"[%s%s]: %s %.2s ap:%d beam:%d")
call pargstr (IMNAME(sh))
call pargstr (IMSEC(sh))
call pargstr (TITLE(sh))
call pargr (IT(sh))
call pargi (AP(sh))
call pargi (BEAM(sh))

# Set GTOOLS labels.
call gt_sets (gt, GTTITLE, Memc[str])
if (UN_LABEL(UN(sh)) != EOS) {
call gt_sets (gt, GTXLABEL, UN_LABEL(UN(sh)))
call gt_sets (gt, GTXUNITS, UN_UNITS(UN(sh)))
} else {
call gt_sets (gt, GTXLABEL, LABEL(sh))
call gt_sets (gt, GTXUNITS, UNITS(sh))
}

call sfree (sp)
end
#@(#) shdr.x 1.1
include
include
include
include "smw.h"
include "units.h"
include "funits.h"
include $$


# SHDR_OPEN -- Open the SHDR spectrum header structure.
# SHDR_GTYPE -- Get the selected spectrum type.
# SHDR_CLOSE -- Close and free the SHDR structure.
# SHDR_COPY -- Make a copy of an SHDR structure.
# SHDR_SYSTEM -- Set or change the WCS system.
# SHDR_UNITS -- Set or change user units.
# SHDR_LW -- Logical to world coordinate transformation.
# SHDR_WL -- World to logical coordinate transformation.
# SHDR_REBIN -- Rebin spectrum to dispersion of reference spectrum.
# SHDR_LINEAR -- Rebin spectrum to linear dispersion.
# SHDR_EXTRACT -- Extract a specific wavelength region.
# SHDR_GI -- Load an integer value from the header.
# SHDR_GR -- Load a real value from the header.


# SHDR_OPEN -- Open SHDR spectrum header structure.
#
# This routine sets header information, WCS transformations, and extracts the
# spectrum from EQUISPEC, MULTISPEC, and NDSPEC format images. The spectrum
# from a 2D/3D format is specified by a logical line and band number.
# Optionally an EQUISPEC or MULTISPEC spectrum may be selected by it's
# aperture number. The access modes are header only or header and data.
# Special checks are made to avoid repeated setting of the header and WCS
# information common to all spectra in an image provided the previously set
# structure is input. Note that the logical to world and world to logical
# transformations require that the MWCS pointer not be closed.

procedure shdr_open (im, smw, index1, index2, ap, mode, sh)

pointer im # IMIO pointer
pointer smw # SMW pointer
int index1 # Image index desired
int index2 # Image index desired
int ap # Aperture number desired
int mode # Access mode
pointer sh # SHDR pointer

int i, j, k, l, n, np, np1, np2, aplow[2], aphigh[2]
real smw_c1tranr(), asumr()
double dval, shdr_lw()
bool newim, streq()
pointer sp, key, str, coeff, mw, ct
pointer smw_sctran(), imgs3r(), un_open(), fun_open()
errchk smw_sctran, imgstr, imgeti, imgetr, smw_gwattrs
errchk un_open, fun_open, fun_ctranr

define data_ 90

begin
call smark (sp)
call salloc (key, SZ_FNAME, TY_CHAR)
call salloc (str, SZ_LINE, TY_CHAR)

# Allocate basic structure or check if the same spectrum is requested.
if (sh == NULL) {
call calloc (sh, LEN_SHDR, TY_STRUCT)
call calloc (SID(sh,1), LEN_SHDRS, TY_CHAR)
newim = true
} else {
call imstats (im, IM_IMAGENAME, Memc[str], SZ_LINE)
newim = !streq (Memc[str], IMNAME(sh))
if (!newim) {
if (LINDEX(sh,1)==index1 && LINDEX(sh,2)==index2) {
if (IS_INDEFI(ap) || AP(sh)==ap) {
np1 = NP1(sh)
np2 = NP2(sh)
np = np2 - np1 + 1
if (CTLW(sh) == NULL || CTWL(sh) == NULL)
goto data_
if (mode == SHHDR) {
do i = 1, SH_NTYPES
call mfree (SPEC(sh,i), TY_REAL)
} else {
switch (SMW_FORMAT(smw)) {
case SMW_ND:
if (mode == SHDATA && SPEC(sh,mode) == NULL)
goto data_
case SMW_ES, SMW_MS:
if (SPEC(sh,mode) == NULL)
goto data_
}
}
call sfree (sp)
return
}
}
}
}

# Set parameters common to an entire image.
if (newim) {
call imstats (im, IM_IMAGENAME, IMNAME(sh), LEN_SHDRS)
IM(sh) = im
MW(sh) = smw

# Get standard parameters.
call shdr_gi (im, "OFLAG", OBJECT, OFLAG(sh))
call shdr_gr (im, "EXPOSURE", INDEF, IT(sh))
call shdr_gr (im, "ITIME", IT(sh), IT(sh))
call shdr_gr (im, "EXPTIME", IT(sh), IT(sh))
call shdr_gr (im, "RA", INDEF, RA(sh))
call shdr_gr (im, "DEC", INDEF, DEC(sh))
call shdr_gr (im, "UT", INDEF, UT(sh))
call shdr_gr (im, "ST", INDEF, ST(sh))
call shdr_gr (im, "HA", INDEF, HA(sh))
call shdr_gr (im, "AIRMASS", INDEF, AM(sh))
call shdr_gi (im, "DC-FLAG", DCNO, DC(sh))
call shdr_gi (im, "EX-FLAG", ECNO, EC(sh))
call shdr_gi (im, "CA-FLAG", FCNO, FC(sh))
iferr (call imgstr (im, "DEREDDEN", RC(sh), LEN_SHDRS))
RC(sh) = EOS

# Flag bad airmass value; i.e. 0.
if (!IS_INDEF (AM(sh)) && AM(sh) < 1.)
AM(sh) = INDEF

# Set the SMW information.
if (SMW_FORMAT(smw) == SMW_MS)
i = 3B
else
i = 2 ** (SMW_PAXIS(smw,1) - 1)
CTLW1(sh) = smw_sctran (smw, "logical", "world", i)
CTWL1(sh) = smw_sctran (smw, "world", "logical", i)

# Set the units.
mw = SMW_MW(smw,0)
i = SMW_PAXIS(smw,1)
iferr (call mw_gwattrs (mw, i, "label", LABEL(sh),LEN_SHDRS))
call strcpy ("", LABEL(sh), LEN_SHDRS)
if (streq (LABEL(sh), "equispec") || streq (LABEL(sh), "multispe"))
call strcpy ("", LABEL(sh), LEN_SHDRS)
iferr (call mw_gwattrs (mw, i, "units", UNITS(sh),LEN_SHDRS))
call strcpy ("", UNITS(sh), LEN_SHDRS)
MWUN(sh) = un_open (UNITS(sh))
call un_copy (MWUN(sh), UN(sh))

iferr (call imgstr (im, "bunit", Memc[str], SZ_LINE))
call strcpy ("", Memc[str], SZ_LINE)
FUNIM(sh) = fun_open (Memc[str])
if (FUN_CLASS(FUNIM(sh)) != FUN_UNKNOWN)
FC(sh) = FCYES

call fun_copy (FUNIM(sh), FUN(sh))
call strcpy (FUN_LABEL(FUN(sh)), FLABEL(sh), LEN_SHDRS)
call strcpy (FUN_UNITS(FUN(sh)), FUNITS(sh), LEN_SHDRS)
}

# Set WCS parameters for spectrum type.
switch (SMW_FORMAT(smw)) {
case SMW_ND:
# Set physical and logical indices.
if (!IS_INDEFI (ap)) {
i = max (1, min (SMW_NSPEC(smw), ap))
j = 1
} else {
i = max (1, index1)
j = max (1, index2)
}
call smw_mw (smw, i, j, mw, k, l)

LINDEX(sh,1) = max (1, min (SMW_LLEN(smw,2), k))
LINDEX(sh,2) = max (1, min (SMW_LLEN(smw,3), l))
PINDEX(sh,1) = LINDEX(sh,1)
PINDEX(sh,2) = LINDEX(sh,2)
APINDEX(sh) = LINDEX(sh,1)

# Set aperture information. Note the use of the logical index.
np1 = 1
call smw_gwattrs (smw, i, j, AP(sh), BEAM(sh), DC(sh),
dval, dval, np2, dval, APLOW(sh,1), APHIGH(sh,1), coeff)

call smw_gapid (smw, i, j, TITLE(sh), LEN_SHDRS)
Memc[SID(sh,1)] = EOS

switch (SMW_LDIM(smw)) {
case 1:
IMSEC(sh) = EOS
case 2:
if (APLOW(sh,1) == APHIGH(sh,1)) {
if (SMW_LAXIS(smw,1) == 1)
call sprintf (IMSEC(sh), LEN_SHDRS, "[*,%d]")
else
call sprintf (IMSEC(sh), LEN_SHDRS, "[%d,*]")
call pargi (nint (APLOW(sh,1)))
} else {
if (SMW_LAXIS(smw,1) == 1)
call sprintf (IMSEC(sh), LEN_SHDRS, "[*,%d:%d]")
else
call sprintf (IMSEC(sh), LEN_SHDRS, "[%d:%d,*]")
call pargi (nint (APLOW(sh,1)))
call pargi (nint (APHIGH(sh,1)))
}
case 3:
if (APLOW(sh,1)==APHIGH(sh,1) && APLOW(sh,2)==APHIGH(sh,2)) {
switch (SMW_LAXIS(smw,1)) {
case 1:
call sprintf (IMSEC(sh), LEN_SHDRS, "[*,%d,%d]")
case 2:
call sprintf (IMSEC(sh), LEN_SHDRS, "[%d,*,%d]")
case 3:
call sprintf (IMSEC(sh), LEN_SHDRS, "[%d,%d,*]")
}
call pargi (nint (APLOW(sh,1)))
call pargi (nint (APLOW(sh,2)))
} else if (APLOW(sh,1) == APHIGH(sh,1)) {
switch (SMW_LAXIS(smw,1)) {
case 1:
call sprintf (IMSEC(sh), LEN_SHDRS, "[*,%d,%d:%d]")
case 2:
call sprintf (IMSEC(sh), LEN_SHDRS, "[%d,*,%d:%d]")
case 3:
call sprintf (IMSEC(sh), LEN_SHDRS, "[%d,%d:%d,*]")
}
call pargi (nint (APLOW(sh,1)))
call pargi (nint (APLOW(sh,2)))
call pargi (nint (APHIGH(sh,2)))
} else if (APLOW(sh,2) == APHIGH(sh,2)) {
switch (SMW_LAXIS(smw,1)) {
case 1:
call sprintf (IMSEC(sh), LEN_SHDRS, "[*,%d:%d,%d]")
case 2:
call sprintf (IMSEC(sh), LEN_SHDRS, "[%d:%d,*,%d]")
case 3:
call sprintf (IMSEC(sh), LEN_SHDRS, "[%d:%d,%d,*]")
}
call pargi (nint (APLOW(sh,1)))
call pargi (nint (APHIGH(sh,1)))
call pargi (nint (APLOW(sh,2)))
} else {
switch (SMW_LAXIS(smw,1)) {
case 1:
call sprintf (IMSEC(sh), LEN_SHDRS, "[*,%d:%d,%d:%d]")
case 2:
call sprintf (IMSEC(sh), LEN_SHDRS, "[%d:%d,*,%d:%d]")
case 3:
call sprintf (IMSEC(sh), LEN_SHDRS, "[%d:%d,%d:%d,*]")
}
call pargi (nint (APLOW(sh,1)))
call pargi (nint (APHIGH(sh,1)))
call pargi (nint (APLOW(sh,2)))
call pargi (nint (APHIGH(sh,2)))
}
}

case SMW_ES, SMW_MS:
# Set the image and aperture indices.
if (SMW_PAXIS(smw,2) != 3) {
PINDEX(sh,1) = max (1, min (SMW_LLEN(smw,2), index1))
PINDEX(sh,2) = max (1, min (SMW_LLEN(smw,3), index2))
LINDEX(sh,1) = PINDEX(sh,1)
LINDEX(sh,2) = PINDEX(sh,2)
APINDEX(sh) = LINDEX(sh,1)
} else {
PINDEX(sh,1) = 1
PINDEX(sh,2) = max (1, min (SMW_LLEN(smw,2), index2))
LINDEX(sh,1) = PINDEX(sh,2)
LINDEX(sh,2) = 1
APINDEX(sh) = 1
}

# If an aperture is specified first try and find it.
# If it is not specified or found then use the physical index.

coeff = NULL
AP(sh) = 0
if (!IS_INDEFI(ap)) {
do i = 1, SMW_NSPEC(smw) {
call smw_gwattrs (smw, i, 1, AP(sh), BEAM(sh), DC(sh),
dval, dval, np2, dval, APLOW(sh,1), APHIGH(sh,1), coeff)
if (AP(sh) == ap && SMW_PAXIS(smw,2) != 3) {
PINDEX(sh,1) = i
LINDEX(sh,1) = i
APINDEX(sh) = i
break
}
}
}
if (AP(sh) != ap)
call smw_gwattrs (smw, APINDEX(sh), 1, AP(sh), BEAM(sh), DC(sh),
dval, dval, np2, dval, APLOW(sh,1), APHIGH(sh,1), coeff)
call mfree (coeff, TY_CHAR)

np1 = 1
if (SMW_PDIM(smw) > 1) {
ct = smw_sctran (smw, "logical", "physical", 2)
PINDEX(sh,1) = nint (smw_c1tranr (ct, real(PINDEX(sh,1))))
call smw_ctfree (ct)
}
if (SMW_PDIM(smw) > 2) {
ct = smw_sctran (smw, "logical", "physical", 4)
PINDEX(sh,2) = nint (smw_c1tranr (ct, real(PINDEX(sh,2))))
call smw_ctfree (ct)
}

call smw_gapid (smw, APINDEX(sh), 1, TITLE(sh), LEN_SHDRS)
call sprintf (Memc[key], SZ_FNAME, "BANDID%d")
call pargi (PINDEX(sh,2))
iferr (call imgstr (im, Memc[key], Memc[SID(sh,1)], LEN_SHDRS))
Memc[SID(sh,1)] = EOS

switch (SMW_LDIM(smw)) {
case 1:
IMSEC(sh) = EOS
case 2:
call sprintf (IMSEC(sh), LEN_SHDRS, "[*,%d]")
call pargi (APINDEX(sh))
case 3:
call sprintf (IMSEC(sh), LEN_SHDRS, "[*,%d,%d]")
call pargi (APINDEX(sh))
call pargi (LINDEX(sh,2))
}
}

# Set NP1 and NP2 in logical coordinates.
i = 2 ** (SMW_PAXIS(smw,1) - 1)
ct = smw_sctran (smw, "physical", "logical", i)
i = max (1, min (int (smw_c1tranr (ct, real (np1))), SMW_LLEN(smw,1)))
j = max (1, min (int (smw_c1tranr (ct, real (np2))), SMW_LLEN(smw,1)))
call smw_ctfree (ct)
np1 = min (i, j)
np2 = max (i, j)
np = np2 - np1 + 1

NP1(sh) = np1
NP2(sh) = np2
SN(sh) = np


data_ # Set the coordinate and data arrays if desired otherwise free them.
CTLW(sh) = CTLW1(sh)
CTWL(sh) = CTWL1(sh)

# Set linear dispersion terms.
W0(sh) = shdr_lw (sh, double(1))
W1(sh) = shdr_lw (sh, double(np))
WP(sh) = (W1(sh) - W0(sh)) / (np - 1)
SN(sh) = np

if (mode == SHHDR) {
do i = 1, SH_NTYPES
call mfree (SPEC(sh,i), TY_REAL)
call sfree (sp)
return
}

# Set WCS array
if (SX(sh) == NULL)
call malloc (SX(sh), np, TY_REAL)
else
call realloc (SX(sh), np, TY_REAL)
do i = 1, np
Memr[SX(sh)+i-1] = shdr_lw (sh, double(i))

# Set spectrum array in most efficient way.
switch (SMW_FORMAT(smw)) {
case SMW_ND:
if (mode == SHDATA || SY(sh) == NULL) {
if (SY(sh) == NULL)
call malloc (SY(sh), np, TY_REAL)
else
call realloc (SY(sh), np, TY_REAL)
call aclrr (Memr[SY(sh)], np)
if (IS_INDEF(APLOW(sh,1)))
aplow[1] = 1
else
aplow[1] = nint (APLOW(sh,1))
if (IS_INDEF(APHIGH(sh,1)))
aphigh[1] = 1
else
aphigh[1] = nint (APHIGH(sh,1))
if (IS_INDEF(APLOW(sh,2)))
aplow[2] = 1
else
aplow[2] = nint (APLOW(sh,2))
if (IS_INDEF(APHIGH(sh,2)))
aphigh[2] = 1
else
aphigh[2] = nint (APHIGH(sh,2))
k = aplow[1]
l = aphigh[1]
n = aphigh[1] - aplow[1] + 1
if (SMW_LAXIS(smw,1) == 1) {
do j = aplow[2], aphigh[2]
do i = aplow[1], aphigh[1]
call aaddr (Memr[imgs3r(im,np1,np2,i,i,j,j)],
Memr[SY(sh)], Memr[SY(sh)], np)
} else if (SMW_LAXIS(smw,1) == 2) {
do j = aplow[2], aphigh[2]
do i = np1, np2
Memr[SY(sh)+i-np1] = Memr[SY(sh)+i-np1] +
asumr (Memr[imgs3r(im,k,l,i,i,j,j)], n)
} else {
do i = np1, np2
do j = aplow[2], aphigh[2]
Memr[SY(sh)+i-np1] = Memr[SY(sh)+i-np1] +
asumr (Memr[imgs3r(im,k,l,j,j,i,i)], n)
}
}
case SMW_ES, SMW_MS:
if (mode == SHDATA || SY(sh) == NULL) {
if (SY(sh) == NULL)
call malloc (SY(sh), np, TY_REAL)
else
call realloc (SY(sh), np, TY_REAL)
i = LINDEX(sh,1)
j = LINDEX(sh,2)
call amovr (Memr[imgs3r(im,np1,np2,i,i,j,j)], Memr[SY(sh)], np)
}

if (mode > SHDATA)
call shdr_gtype (sh, mode)
}

# Guess flux scale if necessary.
if (FC(sh) == FCYES && FUN_CLASS(FUNIM(sh)) == FUN_UNKNOWN) {
if (Memr[SY(sh)+np/2] < 1e-18)
call strcpy ("erg/cm2/s/Hz", Memc[str], SZ_LINE)
else if (Memr[SY(sh)+np/2] < 1e-5)
call strcpy ("erg/cm2/s/A", Memc[str], SZ_LINE)
call fun_close (FUNIM(sh))
FUNIM(sh) = fun_open (Memc[str])
if (FUN_CLASS(FUN(sh)) == FUN_UNKNOWN) {
call fun_copy (FUNIM(sh), FUN(sh))
call strcpy (FUN_LABEL(FUN(sh)), FLABEL(sh), LEN_SHDRS)
call strcpy (FUN_UNITS(FUN(sh)), FUNITS(sh), LEN_SHDRS)
}
}
iferr (call fun_ctranr (FUNIM(sh), FUN(sh), UN(sh), Memr[SX(sh)],
Memr[SY(sh)], Memr[SY(sh)], SN(sh)))
;

call sfree (sp)
end


# SHDR_GTYPE -- Get the selected spectrum type.
# Currently this only works for multispec data.

procedure shdr_gtype (sh, type)

pointer sh #I SHDR pointer
int type #I Spectrum type

int i, j, id, ctowrd(), strdic()
pointer sp, key, str, im, smw, ct, smw_sctran(), imgs3r()
real smw_c1tranr()

begin
im = IM(sh)
smw = MW(sh)

if (SMW_FORMAT(smw) == SMW_ND)
return
if (SMW_PDIM(smw) < 3) {
if (type != SHDATA && type != SHRAW) {
if (SID(sh,type) != NULL)
call mfree (SID(sh,type), TY_CHAR)
if (SPEC(sh,type) != NULL)
call mfree (SPEC(sh,type), TY_REAL)
}
return
}

# Find the band.
call smark (sp)
call salloc (key, SZ_LINE, TY_CHAR)
call salloc (str, SZ_LINE, TY_CHAR)
call salloc (id, SZ_LINE, TY_CHAR)

switch (type) {
case SHDATA:
call strcpy ("|spectrum|data|", Memc[id], SZ_LINE)
case SHRAW:
call strcpy ("|raw|", Memc[id], SZ_LINE)
case SHSIG:
call strcpy ("|sigma|", Memc[id], SZ_LINE)
case SHSKY:
call strcpy ("|background|", Memc[id], SZ_LINE)
case SHCONT:
call strcpy ("|continuum|", Memc[id], SZ_LINE)
}

do i = 1, 4 {
call sprintf (Memc[key], SZ_LINE, "BANDID%d")
call pargi (i)
ifnoerr (call imgstr (im, Memc[key], Memc[str], SZ_LINE)) {
j = 1
if (ctowrd (Memc[str], j, Memc[key], SZ_LINE) == 0)
next
if (strdic (Memc[key], Memc[key], SZ_LINE, Memc[id]) == 0)
next
if (SID(sh,type) == NULL)
call malloc (SID(sh,type), LEN_SHDRS, TY_CHAR)
call strcpy (Memc[str], Memc[SID(sh,type)], LEN_SHDRS)
break
}
}
call sfree (sp)
if (i == 5) {
if (SID(sh,type) != NULL)
call mfree (SID(sh,type), TY_CHAR)
if (SPEC(sh,type) != NULL)
call mfree (SPEC(sh,type), TY_REAL)
return
}

# Map the physical band to logical vector.
ct = smw_sctran (smw, "physical", "logical", 4)
i = nint (smw_c1tranr (ct, real(i)))
call smw_ctfree (ct)
if (SMW_PAXIS(smw,2) != 3) {
if (i > SMW_LLEN(smw,3))
return
j = i
i = LINDEX(sh,1)
} else {
if (i > SMW_LLEN(smw,2))
return
j = 1
}

# Get the spectrum.
if (SPEC(sh,type) == NULL)
call malloc (SPEC(sh,type), SN(sh), TY_REAL)
else
call realloc (SPEC(sh,type), SN(sh), TY_REAL)
call amovr (Memr[imgs3r(im,NP1(sh),NP2(sh),i,i,j,j)],
Memr[SPEC(sh,type)], SN(sh))
end


# SHDR_CLOSE -- Close and free the SHDR structure.

procedure shdr_close (sh)

pointer sh # SHDR structure
int i

begin
if (sh == NULL)
return
do i = 1, SH_NTYPES {
call mfree (SPEC(sh,i), TY_REAL)
call mfree (SID(sh,i), TY_CHAR)
}
call un_close (UN(sh))
call un_close (MWUN(sh))
call fun_close (FUN(sh))
call fun_close (FUNIM(sh))
if (MW(sh) != NULL) {
call smw_ctfree (CTLW1(sh))
call smw_ctfree (CTWL1(sh))
}
call mfree (sh, TY_STRUCT)
end


# SHDR_COPY -- Make a copy of an SHDR structure.
# The image pointer is not copied and the MWCS pointer and transform pointers
# may or may not be copied . The uncopied pointers mean that they will be
# shared by multiple spectrum structures but it also means that when they are
# closed the structures will have invalid pointers. The advantage of not
# copying is that many spectra may come from the same image and the overhead
# of having copies of the IMIO and MWCS pointers can be avoided.

procedure shdr_copy (sh1, sh2, wcs)

pointer sh1 # SHDR structure to copy
pointer sh2 # SHDR structure copy
int wcs # Make copy of wcs?

int i
pointer un, mwun, fun, funim, spec[SH_NTYPES], sid[SH_NTYPES], smw_newcopy()
errchk shdr_system

begin
if (sh2 == NULL) {
call calloc (sh2, LEN_SHDR, TY_STRUCT)
call calloc (SID(sh2,1), LEN_SHDRS, TY_CHAR)
}

un = UN(sh2)
mwun = MWUN(sh2)
fun = FUN(sh2)
funim = FUNIM(sh2)
call amovi (SPEC(sh2,1), spec, SH_NTYPES)
call amovi (SID(sh2,1), sid, SH_NTYPES)
call amovi (Memi[sh1], Memi[sh2], LEN_SHDR)
call amovi (spec, SPEC(sh2,1), SH_NTYPES)
call amovi (sid, SID(sh2,1), SH_NTYPES)
UN(sh2) = un
MWUN(sh2) = mwun
FUN(sh2) = fun
FUNIM(sh2) = funim
call un_copy (UN(sh1), UN(sh2))
call un_copy (MWUN(sh1), MWUN(sh2))
call fun_copy (FUN(sh1), FUN(sh2))
call fun_copy (FUNIM(sh1), FUNIM(sh2))
do i = 1, SH_NTYPES {
if (SPEC(sh1,i) != NULL) {
if (SPEC(sh2,i) == NULL)
call malloc (SPEC(sh2,i), SN(sh1), TY_REAL)
else
call realloc (SPEC(sh2,i), SN(sh1), TY_REAL)
call amovr (Memr[SPEC(sh1,i)], Memr[SPEC(sh2,i)], SN(sh1))
}
}

if (wcs == YES && MW(sh1) != NULL) {
MW(sh2) = smw_newcopy (MW(sh1))
CTLW1(sh2) = NULL
CTWL1(sh2) = NULL
call shdr_system (sh2, "world")
}
end


# SHDR_SYSTEM -- Set or change the WCS system.

procedure shdr_system (sh, system)

pointer sh # SHDR pointer
char system[ARB] # System

int i, sn
bool streq()
double shdr_lw()
pointer smw, mw, smw_sctran(), un_open()
errchk smw_sctran, un_open

begin
smw = MW(sh)
if (smw == NULL)
call error (1, "shdr_system: MWCS not defined")

call smw_ctfree (CTLW1(sh))
call smw_ctfree (CTWL1(sh))

switch (SMW_FORMAT(smw)) {
case SMW_ND, SMW_ES:
i = 2 ** (SMW_PAXIS(smw,1) - 1)
CTLW1(sh) = smw_sctran (smw, "logical", system, i)
CTWL1(sh) = smw_sctran (smw, system, "logical", i)
case SMW_MS:
CTLW1(sh) = smw_sctran (smw, "logical", system, 3B)
CTWL1(sh) = smw_sctran (smw, system, "logical", 3B)
}
CTLW(sh) = CTLW1(sh)
CTWL(sh) = CTWL1(sh)

# Set labels and units
call un_close (MWUN(sh))
if (streq (system, "physical")) {
call strcpy ("Pixel", LABEL(sh), LEN_SHDRS)
call strcpy ("", UNITS(sh), LEN_SHDRS)
MWUN(sh) = un_open (UNITS(sh))
} else {
call smw_mw (smw, 1, 1, mw, i, i)
iferr (call mw_gwattrs (mw, SMW_PAXIS(smw,1), "label", LABEL(sh),
LEN_SHDRS))
call strcpy ("", LABEL(sh), LEN_SHDRS)
if (streq (LABEL(sh), "equispec") || streq (LABEL(sh), "multispe"))
call strcpy ("", LABEL(sh), LEN_SHDRS)
iferr (call mw_gwattrs (mw, SMW_PAXIS(smw,1), "units", UNITS(sh),
LEN_SHDRS))
call strcpy ("", UNITS(sh), LEN_SHDRS)
MWUN(sh) = un_open (UNITS(sh))
call strcpy (UN_LABEL(UN(sh)), LABEL(sh), LEN_SHDRS)
call strcpy (UN_UNITS(UN(sh)), UNITS(sh), LEN_SHDRS)
}

sn = SN(sh)
W0(sh) = shdr_lw (sh, double(1))
W1(sh) = shdr_lw (sh, double(sn))
WP(sh) = (W1(sh) - W0(sh)) / (sn - 1)
if (SX(sh) != NULL)
do i = 1, sn
Memr[SX(sh)+i-1] = shdr_lw (sh, double(i))
end


# SHDR_UNITS -- Set or change the WCS system.
# This changes W0, W1, WP, and SX.

procedure shdr_units (sh, units)

pointer sh # SHDR pointer
char units[ARB] # Units

int i, sn
bool streq()
double shdr_lw()
pointer str, un, un_open()
errchk un_open

begin
# Check for unknown units.
if (streq (units, "display")) {
call malloc (str, SZ_LINE, TY_CHAR)
iferr (call mw_gwattrs (SMW_MW(MW(sh),0), SMW_PAXIS(MW(sh),1),
"units_display", Memc[str], SZ_LINE)) {
un = NULL
call un_copy (MWUN(sh), un)
} else
un = un_open (Memc[str])
call mfree (str, TY_CHAR)
} else if (streq (units, "default")) {
un = NULL
call un_copy (MWUN(sh), un)
} else
un = un_open (units)
if (UN_CLASS(un) == UN_UNKNOWN || UN_CLASS(MWUN(sh)) == UN_UNKNOWN) {
call un_close (un)
call error (1, "Cannot convert to specified units")
}

# Update the coordinates.
call un_close (UN(sh))
UN(sh) = un

call strcpy (UN_LABEL(UN(sh)), LABEL(sh), LEN_SHDRS)
call strcpy (UN_UNITS(UN(sh)), UNITS(sh), LEN_SHDRS)

sn = SN(sh)
W0(sh) = shdr_lw (sh, double(1))
W1(sh) = shdr_lw (sh, double(sn))
WP(sh) = (W1(sh) - W0(sh)) / (sn - 1)
if (SX(sh) != NULL)
do i = 1, sn
Memr[SX(sh)+i-1] = shdr_lw (sh, double(i))
end


# SHDR_LW -- Logical to world coordinate transformation.
# The transformation pointer is generally NULL only after SHDR_LINEAR

double procedure shdr_lw (sh, l)

pointer sh # SHDR pointer
double l # Logical coordinate
double w # World coordinate

double l0, l1, l2, w1, smw_c1trand()

begin
l0 = l + NP1(sh) - 1
if (CTLW(sh) != NULL) {
switch (SMW_FORMAT(MW(sh))) {
case SMW_ND, SMW_ES:
w = smw_c1trand (CTLW(sh), l0)
if (DC(sh) == DCLOG)
if (SMW_CTTYPE(CTLW(sh)) == SMW_LW)
w = 10. ** max (-20D0, min (20D0, w))
case SMW_MS:
call smw_c2trand (CTLW(sh), l0, double (APINDEX(sh)), w, w1)
}
} else {
switch (DC(sh)) {
case DCLOG:
w = W0(sh) * 10. ** (log10(W1(sh)/W0(sh)) * (l0-1) / (SN(sh)-1))
case DCFUNC:
w = W0(sh)
call smw_c2trand (CTWL1(sh), w, double (AP(sh)), l1, w1)
w = W1(sh)
call smw_c2trand (CTWL1(sh), w, double (AP(sh)), l2, w1)
if (SN(sh) > 1)
l1 = (l2 - l1) / (SN(sh) - 1) * (l0 - 1) + l1
else
l1 = l0 - 1 + l1
call smw_c2trand (CTLW1(sh), l1, double (APINDEX(sh)), w, w1)
default:
w = W0(sh) + (l0 - 1) * WP(sh)
}
}

iferr (call un_ctrand (MWUN(sh), UN(sh), w, w, 1))
;
return (w)
end


# SHDR_WL -- World to logical coordinate transformation.
# The transformation pointer is generally NULL only after SHDR_LINEAR

double procedure shdr_wl (sh, w)

pointer sh # SHDR pointer
double w # World coordinate
double l # Logical coordinate

double w1, l1, l2, smw_c1trand()

begin
iferr (call un_ctrand (UN(sh), MWUN(sh), w, w1, 1))
w1 = w

if (CTWL(sh) != NULL) {
switch (SMW_FORMAT(MW(sh))) {
case SMW_ND, SMW_ES:
if (DC(sh) == DCLOG)
if (SMW_CTTYPE(CTWL(sh)) == SMW_WL)
w1 = log10 (w1)
l = smw_c1trand (CTWL(sh), w1)
case SMW_MS:
call smw_c2trand (CTWL(sh), w1, double (AP(sh)),l,l1)
}
} else {
switch (DC(sh)) {
case DCLOG:
l = log10(w1/W0(sh)) / log10(W1(sh)/W0(sh)) * (SN(sh)-1) + 1
case DCFUNC:
call smw_c2trand (CTWL1(sh), w1, double (AP(sh)), l, l1)

w1 = W0(sh)
call smw_c2trand (CTWL1(sh), w1, double (AP(sh)), l1, w1)
w1 = W1(sh)
call smw_c2trand (CTWL1(sh), w1, double (AP(sh)), l2, w1)
if (l1 != l2)
l = (SN(sh) - 1) / (l2 - l1) * (l - l1) + 1
else
l = l - l1 + 1
default:
l = (w1 - W0(sh)) / WP(sh) + 1
}
}

return (l-NP1(sh)+1)
end


# SHDR_REBIN -- Rebin spectrum to dispersion of reference spectrum.
# The interpolation function is set by ONEDINTERP.

procedure shdr_rebin (sh, shref)

pointer sh # Spectrum to be rebinned
pointer shref # Reference spectrum

char interp[10]
int i, j, type, ia, ib, n, clgwrd()
real a, b, sum, asieval(), asigrl()
double x, w, xmin, xmax, shdr_lw(), shdr_wl()
pointer unsave, asi, spec
bool fp_equalr()

begin
# Check if rebinning is needed
if (DC(sh) == DC(shref) && DC(sh) != DCFUNC &&
fp_equalr (W0(sh), W0(shref)) && fp_equalr(WP(sh), WP(shref)) &&
SN(sh) == SN(shref))
return

# Do everything in units of MWCS.
unsave = UN(sh)
UN(SH) = MWUN(sh)

call asiinit (asi, clgwrd ("interp", interp, 10, II_FUNCTIONS))
do type = 1, SH_NTYPES {
if (SPEC(sh,type) == NULL)
next

# Fit the interpolation function to the spectrum.
# Extend the interpolation by one pixel at each end.

n = SN(sh)
call malloc (spec, n+2, TY_REAL)
call amovr (Memr[SPEC(sh,type)], Memr[spec+1], n)
Memr[spec] = Memr[SPEC(sh,type)]
Memr[spec+n+1] = Memr[SPEC(sh,type)+n-1]
call asifit (asi, Memr[spec], n+2)
call mfree (spec, TY_REAL)

xmin = 0.5
xmax = n + 0.5

# Reallocate data array
if (n != SN(shref)) {
n = SN(shref)
call realloc (SPEC(sh,type), n, TY_REAL)
call aclrr (Memr[SPEC(sh,type)], n)
}
spec = SPEC(sh,type)

# Compute the average flux in each output pixel.

x = 0.5
w = shdr_lw (shref, x)
x = shdr_wl (sh, w)
b = max (xmin, min (xmax, x)) + 1
do i = 1, n {
x = i + 0.5
w = shdr_lw (shref, x)
x = shdr_wl (sh, w)
a = b
b = max (xmin, min (xmax, x)) + 1
if (a <= b) {
ia = nint (a + 0.5)
ib = nint (b - 0.5)
if (abs (a+0.5-ia) < .001 && abs (b-0.5-ib) < .001) {
sum = 0.
do j = ia, ib
sum = sum + asieval (asi, real(j))
if (ib - ia > 0)
sum = sum / (ib - ia)
} else {
sum = asigrl (asi, a, b)
if (b - a > 0.)
sum = sum / (b - a)
}
} else {
ib = nint (b + 0.5)
ia = nint (a - 0.5)
if (abs (a-0.5-ia) < .001 && abs (b+0.5-ib) < .001) {
sum = 0.
do j = ib, ia
sum = sum + asieval (asi, real(j))
if (ia - ib > 0)
sum = sum / (ia - ib)
} else {
sum = asigrl (asi, b, a)
if (a - b > 0.)
sum = sum / (a - b)
}
}

Memr[spec] = sum
spec = spec + 1
}
}
call asifree (asi)

# Set the rest of the header. The coordinate transformations are
# canceled to indicate they are not valid for the data. They
# are not freed because the same pointer may be used in other
# spectra from the same image.

if (SN(sh) != n)
call realloc (SX(sh), n, TY_REAL)
call amovr (Memr[SX(shref)], Memr[SX(sh)], n)
DC(sh) = DC(shref)
W0(sh) = W0(shref)
W1(sh) = W1(shref)
WP(sh) = WP(shref)
SN(sh) = SN(shref)

CTLW(sh) = NULL
CTWL(sh) = NULL

# Restore original units
UN(sh) = unsave
iferr (call un_ctranr (MWUN(sh), UN(sh), Memr[SX(sh)], Memr[SX(sh)],
SN(sh)))
;
end


# SHDR_LINEAR -- Rebin spectrum to linear dispersion.
# The interpolation function is set by ONEDINTERP

procedure shdr_linear (sh, w0, w1, sn, dc)

pointer sh # Spectrum to be rebinned
real w0 # Wavelength of first logical pixel
real w1 # Wavelength of last logical pixel
int sn # Number of pixels
int dc # Dispersion type (DCLINEAR | DCLOG)

char interp[10]
int i, j, type, ia, ib, n, clgwrd()
real a, b, sum, asieval(), asigrl()
double x, w, w0l, wp, xmin, xmax, shdr_wl()
pointer unsave, asi, spec
bool fp_equalr()

begin
# Check if rebinning is needed
if (DC(sh) == dc && fp_equalr (W0(sh), w0) &&
fp_equalr (W1(sh), w1) && SN(sh) == sn)
return

# Do everything in units of MWCS.
unsave = UN(sh)
UN(SH) = MWUN(sh)

call asiinit (asi, clgwrd ("interp", interp, 10, II_FUNCTIONS))
do type = 1, SH_NTYPES {
if (SPEC(sh,type) == NULL)
next

# Fit the interpolation function to the spectrum.
# Extend the interpolation by one pixel at each end.

n = SN(sh)
call malloc (spec, n+2, TY_REAL)
call amovr (Memr[SPEC(sh,type)], Memr[spec+1], n)
Memr[spec] = Memr[SPEC(sh,type)]
Memr[spec+n+1] = Memr[SPEC(sh,type)+n-1]
call asifit (asi, Memr[spec], n+2)
call mfree (spec, TY_REAL)

xmin = 0.5
xmax = n + 0.5

# Reallocate spectrum data array
if (n != sn) {
n = sn
call realloc (SPEC(sh,type), n, TY_REAL)
}
spec = SPEC(sh,type)

# Integrate across pixels using ASIGRL.

x = 0.5
if (dc == DCLOG) {
w0l = log10 (w0)
wp = (log10 (w1) - log10(w0)) / (n - 1)
w = 10. ** (w0l+(x-1)*wp)
} else {
wp = (w1 - w0) / (n - 1)
w = w0 + (x - 1) * wp
}
x = shdr_wl (sh, w)
b = max (xmin, min (xmax, x)) + 1
do i = 1, n {
x = i + 0.5
if (dc == DCLOG)
w = 10. ** (w0l + (x - 1) * wp)
else
w = w0 + (x - 1) * wp
x = shdr_wl (sh, w)
a = b
b = max (xmin, min (xmax, x)) + 1
if (a <= b) {
ia = nint (a + 0.5)
ib = nint (b - 0.5)
if (abs (a+0.5-ia) < .001 && abs (b-0.5-ib) < .001) {
sum = 0.
do j = ia, ib
sum = sum + asieval (asi, real(j))
if (ib - ia > 0)
sum = sum / (ib - ia)
} else {
sum = asigrl (asi, a, b)
if (b - a > 0.)
sum = sum / (b - a)
}
} else {
ib = nint (b + 0.5)
ia = nint (a - 0.5)
if (abs (a-0.5-ia) < .001 && abs (b+0.5-ib) < .001) {
sum = 0.
do j = ib, ia
sum = sum + asieval (asi, real(j))
if (ia - ib > 0)
sum = sum / (ia - ib)
} else {
sum = asigrl (asi, b, a)
if (a - b > 0.)
sum = sum / (a - b)
}
}
Memr[spec] = sum
spec = spec + 1
}
}
call asifree (asi)

# Set the rest of the header. The coordinate transformations are
# canceled to indicate they are not valid for the data. They
# are not freed because the same pointer may be used in other
# spectra from the same image.

if (SN(sh) != n)
call realloc (SX(sh), n, TY_REAL)
do i = 0, n-1 {
if (dc == DCLOG)
w = 10. ** (w0l + i * wp)
else
w = w0 + i * wp
Memr[SX(sh)+i] = w
}
W0(sh) = w0
W1(sh) = w1
WP(sh) = (w1 - w0) / (sn - 1)
SN(sh) = sn
NP1(sh) = 1
NP2(sh) = sn
DC(sh) = dc

CTLW(sh) = NULL
CTWL(sh) = NULL

# Restore original units
UN(sh) = unsave
iferr (call un_ctranr (MWUN(sh), UN(sh), Memr[SX(sh)], Memr[SX(sh)],
sn))
;
end


# SHDR_EXTRACT -- Extract a specific wavelength region.

procedure shdr_extract (sh, w1, w2, rebin)

pointer sh # SHDR structure
real w1 # Starting wavelength
real w2 # Ending wavelength
bool rebin # Rebin wavelength region?

int i, j, i1, i2, n
double l1, l2
pointer buf
bool fp_equald()
double shdr_wl(), shdr_lw()
errchk shdr_linear, shdr_lw, shdr_wl

begin
l1 = shdr_wl (sh, double (w1))
l2 = shdr_wl (sh, double (w2))
if (fp_equald(l1,l2) || max(l1,l2) < 1 || min (l1,l2) > SN(sh))
call error (1, "No pixels to extract")
l1 = max (1D0, min (double (SN(sh)), l1))
l2 = max (1D0, min (double (SN(sh)), l2))
i1 = nint (l1)
i2 = nint (l2)
n = abs (i2 - i1) + 1

if (rebin) {
l1 = shdr_lw (sh, l1)
l2 = shdr_lw (sh, l2)
if (DC(sh) == DCFUNC)
call shdr_linear (sh, real (l1), real (l2), n, DCLINEAR)
else
call shdr_linear (sh, real (l1), real (l2), n, DC(sh))
} else {
if (i1 == 1 && i2 == SN(sh))
return

if (i1 <= i2) {
do j = 1, SH_NTYPES
if (SPEC(sh,j) != NULL)
call amovr (Memr[SPEC(sh,j)+i1-1], Memr[SPEC(sh,j)], n)
} else {
call malloc (buf, n, TY_REAL)
do j = 1, SH_NTYPES {
do i = i1, i2, -1
Memr[buf+i1-i] = Memr[SPEC(sh,j)+i-1]
call amovr (Memr[buf], Memr[SPEC(sh,j)], n)
}
call mfree (buf, TY_REAL)
}
W0(sh) = Memr[SX(sh)]
W1(sh) = Memr[SX(sh)+n-1]
SN(sh) = n
NP1(sh) = 1
NP2(sh) = n
if (n > 1)
WP(sh) = (W1(sh) - W0(sh)) / (SN(sh) - 1)
CTLW(sh) = NULL
CTWL(sh) = NULL
}
end


# SHDR_GI -- Load an integer value from the header.

procedure shdr_gi (im, field, default, ival)

pointer im
char field[ARB]
int default
int ival

int dummy, imaccf(), imgeti()

begin
ival = default
if (imaccf (im, field) == YES) {
iferr (dummy = imgeti (im, field))
call erract (EA_WARN)
else
ival = dummy
}
end


# SHDR_GR -- Load a real value from the header.

procedure shdr_gr (im, field, default, rval)

pointer im
char field[ARB]
real default
real rval

int imaccf()
real dummy, imgetr()

begin
rval = default
if (imaccf (im, field) == YES) {
iferr (dummy = imgetr (im, field))
call erract (EA_WARN)
else
rval = dummy
}
end
#@(#) smwclose.x 1.1
include "smw.h"


# SMW_CLOSE -- Close the SMW data structure.
# This includes closing the MWCS pointers.

procedure smw_close (smw)

pointer smw # SMW pointer

int i
pointer apids

begin
if (smw == NULL)
return

switch (SMW_FORMAT(smw)) {
case SMW_ND:
call mfree (SMW_APID(smw), TY_CHAR)
call mw_close (SMW_MW(smw,0))
case SMW_ES:
call mfree (SMW_APS(smw), TY_INT)
call mfree (SMW_BEAMS(smw), TY_INT)
call mfree (SMW_APLOW(smw), TY_REAL)
call mfree (SMW_APHIGH(smw), TY_REAL)
call mfree (SMW_APID(smw), TY_CHAR)
apids = SMW_APIDS(smw) - 1
do i = 1, SMW_NSPEC(smw)
call mfree (Memi[apids+i], TY_CHAR)
call mfree (SMW_APIDS(smw), TY_POINTER)
call mw_close (SMW_MW(smw,0))
case SMW_MS:
call mfree (SMW_APS(smw), TY_INT)
call mfree (SMW_BEAMS(smw), TY_INT)
call mfree (SMW_APLOW(smw), TY_REAL)
call mfree (SMW_APHIGH(smw), TY_REAL)
call mfree (SMW_APID(smw), TY_CHAR)
apids = SMW_APIDS(smw) - 1
do i = 1, SMW_NSPEC(smw)
call mfree (Memi[apids+i], TY_CHAR)
do i = 0, SMW_NMW(smw)-1 {
if (SMW_MW(smw,i) != NULL)
call mw_close (SMW_MW(smw,i))
}
}
call mfree (smw, TY_STRUCT)
end
#@(#) smwct.x 1.1
include "smw.h"


# SMW_CT -- Get MCWS CT pointer for the specified physical line.

pointer procedure smw_ct (sct, line)

pointer sct #I SMW pointer
int line #I Physical line

begin
if (SMW_NCT(sct) == 1)
return (SMW_CT(sct,0))

if (line < 1 || line > SMW_NSPEC(SMW_SMW(sct)))
call error (1, "smw_ct: aperture not found")

return (SMW_CT(sct,(line-1)/SMW_NSPLIT))
end
#@(#) smwctfree.x 1.1
include "smw.h"


# SMW_CTFREE -- Free a spectral SMW coordinate transform pointer.

procedure smw_ctfree (ct)

pointer ct # SMW CT pointer
int i

begin
if (ct == NULL)
return

do i = 0, SMW_NCT(ct)-1
call mw_ctfree (SMW_CT(ct,i))
call mw_ctfree (SMW_CTL(ct))
call mfree (ct, TY_STRUCT)
end
#@(#) smwctran.x 1.2
include "smw.h"


# Evaluate SMW coordinate transform. These procedures simple call the
# MWCS procedures unless the WCS is a split MULTISPEC WCS. In that
# case the appropriate piece needs to be determined and the physical
# line numbers manipulated.
#
# SMW_CTRANR -- N dimensional real coordinate transformation.
# SMW_CTRAND -- N dimensional double coordinate transformation.
# SMW_C1TRANR -- One dimensional real coordinate transformation.
# SMW_C1TRAND -- One dimensional double coordinate transformation.
# SMW_C2TRANR -- Two dimensional real coordinate transformation.
# SMW_C2TRAND -- Two dimensional double coordinate transformation.


# SMW_CTRANR -- N dimensional real coordinate transformation.

procedure smw_ctranr (ct, p1, p2, ndim)

pointer ct #I SMW CT pointer
real p1[ndim] #I Input coordinate
real p2[ndim] #O Output coordinate
int ndim #I Dimensionality

errchk mw_ctranr

begin
if (SMW_NCT(ct) == 1) {
call mw_ctranr (SMW_CT(ct,0), p1, p2, ndim)
return
}

call error (1, "SMW_CTRANR: Wrong WCS type")
end


# SMW_CTRAND -- N dimensional double coordinate transformation.

procedure smw_ctrand (ct, p1, p2, ndim)

pointer ct #I SMW CT pointer
double p1[ndim] #I Input coordinate
double p2[ndim] #O Output coordinate
int ndim #I Dimensionality

errchk mw_ctrand

begin
if (SMW_NCT(ct) == 1) {
call mw_ctrand (SMW_CT(ct,0), p1, p2, ndim)
return
}

call error (1, "SMW_CTRAND: Wrong WCS type")
end


# SMW_C1TRANR -- One dimensional real coordinate transformation.

real procedure smw_c1tranr (ct, x)

pointer ct #I SMW CT pointer
real x #I Input coordinate

real mw_c1tranr()
errchk mw_c1tranr

begin
if (SMW_NCT(ct) == 1)
return (mw_c1tranr (SMW_CT(ct,0), x))

call error (1, "SMW_C1TRANR: Wrong WCS type")
end


# SMW_C1TRAND -- One dimensional double coordinate transformation.

double procedure smw_c1trand (ct, x)

pointer ct #I SMW CT pointer
double x #I Input coordinate

double mw_c1trand()
errchk mw_c1trand

begin
if (SMW_NCT(ct) == 1)
return (mw_c1trand (SMW_CT(ct,0), x))

call error (1, "SMW_C1TRAND: Wrong WCS type")
end


# SMW_C2TRANR -- Two dimensional real coordinate transformation.

procedure smw_c2tranr (ct, x1, y1, x2, y2)

pointer ct #I SMW CT pointer
real x1, y1 #I Input coordinates
real x2, y2 #O Output coordinates

int i, j, naps
real xp, yp
pointer aps, smw_ct()
errchk smw_ct, mw_c2tranr

begin
if (SMW_NCT(ct) == 1) {
# EQUISPEC requires mapping physical lines to apertures.
switch (SMW_FORMAT(SMW_SMW(ct))) {
case SMW_ES:
aps = SMW_APS(SMW_SMW(ct))
naps = SMW_NSPEC(SMW_SMW(ct))
switch (SMW_CTTYPE(ct)) {
case SMW_LW, SMW_PW:
call mw_c2tranr (SMW_CT(ct,0), x1, y1, x2, y2)
i = max (1, min (naps, nint (y2)))
y2 = Memi[aps+i-1]
case SMW_WL, SMW_WP:
j = nint (y1)
yp = 1
do i = 1, naps
if (j == Memi[aps+i-1]) {
yp = i
break
}
call mw_c2tranr (SMW_CT(ct,0), x1, yp, x2, y2)
default:
call mw_c2tranr (SMW_CT(ct,0), x1, y1, x2, y2)
}
default:
call mw_c2tranr (SMW_CT(ct,0), x1, y1, x2, y2)
}
return
}

# If we get here then we are dealing with a split MULTISPEC WCS.
# Depending on the systems being transformed there may need to
# transformation made on the physical coordinate system.

switch (SMW_CTTYPE(ct)) {
case SMW_LW:
call mw_c2tranr (SMW_CTL(ct), x1, y1, xp, yp)
i = nint (yp)
yp = mod (i-1, SMW_NSPLIT) + 1
call mw_c2tranr (smw_ct(ct,i), xp, yp, x2, y2)
case SMW_PW:
i = nint (y1)
yp = mod (i-1, SMW_NSPLIT) + 1
call mw_c2tranr (smw_ct(ct,i), x1, yp, x2, y2)
case SMW_WL:
aps = SMW_APS(SMW_SMW(ct))
naps = SMW_NSPEC(SMW_SMW(ct))
j = nint (y1)
do i = 1, naps {
if (Memi[aps+i-1] == j) {
call mw_c2tranr (smw_ct(ct,i), x1, y1, xp, yp)
yp = i
call mw_c2tranr (SMW_CTL(ct), xp, yp, x2, y2)
return
}
}
call error (1, "Aperture not found")
case SMW_WP:
aps = SMW_APS(SMW_SMW(ct))
naps = SMW_NSPEC(SMW_SMW(ct))
j = nint (y1)
do i = 1, naps {
if (Memi[aps+i-1] == j) {
call mw_c2tranr (smw_ct(ct,i), x1, y1, x2, y2)
y2 = i
return
}
}
call error (1, "Aperture not found")
default:
x2 = x1
y2 = y1
}
end


# SMW_C2TRAND -- Two dimensional double coordinate transformation

procedure smw_c2trand (ct, x1, y1, x2, y2)

pointer ct #I SMW CT pointer
double x1, y1 #I Input coordinates
double x2, y2 #O Output coordinates

int i, j, naps
double xp, yp
pointer aps, smw_ct()
errchk smw_ct, mw_c2trand

begin
if (SMW_NCT(ct) == 1) {
# EQUISPEC requires mapping physical lines to apertures.
switch (SMW_FORMAT(SMW_SMW(ct))) {
case SMW_ES:
aps = SMW_APS(SMW_SMW(ct))
naps = SMW_NSPEC(SMW_SMW(ct))
switch (SMW_CTTYPE(ct)) {
case SMW_LW, SMW_PW:
call mw_c2trand (SMW_CT(ct,0), x1, y1, x2, y2)
i = max (1, min (naps, nint (y2)))
y2 = Memi[aps+i-1]
case SMW_WL, SMW_WP:
j = nint (y1)
yp = 1
do i = 1, naps
if (j == Memi[aps+i-1]) {
yp = i
break
}
call mw_c2trand (SMW_CT(ct,0), x1, yp, x2, y2)
default:
call mw_c2trand (SMW_CT(ct,0), x1, y1, x2, y2)
}
default:
call mw_c2trand (SMW_CT(ct,0), x1, y1, x2, y2)
}
return
}

# If we get here then we are dealing with a split MULTISPEC WCS.
# Depending on the systems being transformed there may need to
# transformation made on the physical coordinate system.

switch (SMW_CTTYPE(ct)) {
case SMW_LW:
call mw_c2trand (SMW_CTL(ct), x1, y1, xp, yp)
i = nint (yp)
yp = mod (i-1, SMW_NSPLIT) + 1
call mw_c2trand (smw_ct(ct,i), xp, yp, x2, y2)
case SMW_PW:
i = nint (y1)
yp = mod (i-1, SMW_NSPLIT) + 1
call mw_c2trand (smw_ct(ct,i), x1, yp, x2, y2)
case SMW_WL:
aps = SMW_APS(SMW_SMW(ct))
naps = SMW_NSPEC(SMW_SMW(ct))
j = nint (y1)
do i = 1, naps {
if (Memi[aps+i-1] == j) {
call mw_c2trand (smw_ct(ct,i), x1, y1, xp, yp)
yp = i
call mw_c2trand (SMW_CTL(ct), xp, yp, x2, y2)
return
}
}
call error (1, "Aperture not found")
case SMW_WP:
aps = SMW_APS(SMW_SMW(ct))
naps = SMW_NSPEC(SMW_SMW(ct))
j = nint (y1)
do i = 1, naps {
if (Memi[aps+i-1] == j) {
call mw_c2trand (smw_ct(ct,i), x1, y1, x2, y2)
y2 = i
return
}
}
call error (1, "Aperture not found")
default:
x2 = x1
y2 = y1
}
end
#@(#) smwdaxis.x 1.1
include "smw.h"

# SMW_DAXIS -- Set physical dispersion axis and summing factors.
# A default value of zero for the dispersion axis will cause the dispersion
# axis to be sought in the image header and, if not found, from the CL
# "dispaxis" parameter. A default value of zero for the summing factors will
# cause them to be queried from the CL "nsum" parameter. A default value of
# INDEFI in either parameter will leave the current default unchanged.
#
# When this procedure is called with an SMW and IMIO pointer the SMW
# pointer is updated to desired default dispersion axis and summing
# parameters.

procedure smw_daxis (smw, im, daxisp, nsum1, nsum2)

pointer smw #I SMW pointer
pointer im #I IMIO pointer
int daxisp #I Default dispersion axis
int nsum1, nsum2 #I Default summing factors

int i, da, ns[2], imgeti(), clgeti(), clscan(), nscan()
data da/0/, ns/0,0/
errchk clgeti

begin
# Set defaults.
# A value of 0 will use the image DISPAXIS or query the CL and
# a value of INDEFI will leave the current default unchanged.

if (!IS_INDEFI (daxisp))
da = daxisp
if (!IS_INDEFI (nsum1))
ns[1] = nsum1
if (!IS_INDEFI (nsum2))
ns[2] = nsum2

if (smw == NULL)
return

# This procedure is specific to the NDSPEC format.
if (SMW_FORMAT(smw) != SMW_ND)
return

# Set dispersion axis.
if (da == 0) {
if (im == NULL)
SMW_PAXIS(smw,1) = clgeti ("dispaxis")
else {
iferr (SMW_PAXIS(smw,1) = imgeti (im, "DISPAXIS"))
SMW_PAXIS(smw,1) = clgeti ("dispaxis")
}
} else
SMW_PAXIS(smw,1) = da

# Set summing parameters.
if (ns[1] == 0 || ns[2] == 0) {
if (clscan("nsum") == EOF)
call error (1, "smw_daxis: Error in 'nsum' parameter")
call gargi (i)
if (ns[1] == 0) {
if (nscan() == 1)
SMW_NSUM(smw,1) = max (1, i)
else
call error (1, "smw_daxis: Error in 'nsum' parameter")
} else
SMW_NSUM(smw,1) = ns[1]
call gargi (i)
if (ns[2] == 0) {
if (nscan() == 2)
SMW_NSUM(smw,2) = max (1, i)
else
SMW_NSUM(smw,2) = SMW_NSUM(smw,1)
} else
SMW_NSUM(smw,2) = ns[2]
} else {
SMW_NSUM(smw,1) = ns[1]
SMW_NSUM(smw,2) = ns[2]
}
end
#@(#) smwequispec.x 1.1
include
include
include "smw.h"


# SMW_EQUISPEC -- Setup the EQUISPEC SMW parameters.
# The aperture information is in the APNUM and APID keywords and the
# WCS information is in the linear MWCS.

procedure smw_equispec (im, smw)

pointer im #I IMIO pointer
pointer smw #U MWCS pointer input SMW pointer output

int i, j, k, nchar, ap, beam, dtype, nw
double w1, dw, z
real aplow[2], aphigh[2], mw_c1tranr()
pointer sp, key, val, wterm, mw, ct, mw_sctran()
int imgeti(), mw_stati(), ctoi(), ctor()
errchk imgstr, mw_gwtermd, mw_sctran
errchk smw_open, smw_saxes, smw_swattrs, smw_sapid

begin
call smark (sp)
call salloc (key, SZ_FNAME, TY_CHAR)
call salloc (val, SZ_LINE, TY_CHAR)

# Determine the dispersion parameters
mw = smw
i = mw_stati (mw, MW_NDIM)
call salloc (wterm, 2*i+i*i, TY_DOUBLE)
call mw_gwtermd (mw, Memd[wterm], Memd[wterm+i], Memd[wterm+2*i], i)
w1 = Memd[wterm+i]
dw = Memd[wterm+2*i]

# Determine the number of physical pixels.
ct = mw_sctran (mw, "logical", "physical", 1)
nw = max (mw_c1tranr (ct, 1.), mw_c1tranr (ct, real(IM_LEN(im,1))))
call mw_ctfree (ct)

# Determine the dispersion type.
iferr (dtype = imgeti (im, "DC-FLAG"))
dtype = DCNO
if (dtype==DCLOG) {
if (abs(w1)>20. || abs(w1+(nw-1)*dw)>20.)
dtype = DCLINEAR
else {
w1 = 10D0 ** w1
dw = w1 * (10D0 ** ((nw-1)*dw) - 1) / (nw - 1)
}
}

# Set the SMW data structure.
call smw_open (smw, NULL, im)
do i = 1, SMW_NSPEC(smw) {
call smw_mw (smw, i, 1, mw, j, k)
call sprintf (Memc[key], SZ_FNAME, "APNUM%d")
call pargi (j)
call imgstr (im, Memc[key], Memc[val], SZ_LINE)
k = 1
nchar = ctoi (Memc[val], k, ap)
nchar = ctoi (Memc[val], k, beam)
if (ctor (Memc[val], k, aplow[1]) == 0)
aplow[1] = INDEF
if (ctor (Memc[val], k, aphigh[1]) == 0)
aphigh[1] = INDEF
if (ctor (Memc[val], k, aplow[2]) == 0)
aplow[2] = INDEF
if (ctor (Memc[val], k, aphigh[2]) == 0)
aphigh[2] = INDEF
z = 0.

call smw_swattrs (smw, i, 1, ap, beam, dtype, w1, dw, nw, z,
aplow, aphigh, "")

call sprintf (Memc[key], SZ_FNAME, "APID%d")
call pargi (j)
ifnoerr (call imgstr (im, Memc[key], Memc[val], SZ_LINE))
call smw_sapid (smw, i, 1, Memc[val])
}

call sfree (sp)
end
#@(#) smwesms.x 1.1
include
include "smw.h"


# SMW_ESMS -- Convert EQUISPEC WCS into MULTISPEC WCS.
# This is called by SMW_SWATTRS when the equal linear coordinate system
# requirement of the EQUISPEC WCS is violated.

procedure smw_esms (smw)

pointer smw #U SMW pointer

int i, j, k, pdim1, pdim2, ap, beam, dtype, nw, axes[2]
double w1, dw, z
real aplow, aphigh
pointer sp, key, str, lterm, mw, mw1, mw2, apid, mw_open()
data axes/1,2/

begin
call smark (sp)
call salloc (key, SZ_FNAME, TY_CHAR)
call salloc (str, SZ_LINE, TY_CHAR)
call salloc (lterm, 12, TY_DOUBLE)

# Set the basic MWCS
mw1 = SMW_MW(smw,0)
pdim1 = SMW_PDIM(smw)
pdim2 = max (2, pdim1)
mw2 = mw_open (NULL, pdim2)
call mw_newsystem (mw2, "multispec", pdim2)
call mw_swtype (mw2, axes, 2, "multispec", "")
if (pdim2 > 2)
call mw_swtype (mw2, 3, 1, "linear", "")
call mw_gltermd (mw1, Memd[lterm+pdim2], Memd[lterm], pdim1)
if (pdim1 == 1) {
Memd[lterm+1] = 0.
Memd[lterm+3] = 0.
Memd[lterm+4] = 0.
Memd[lterm+5] = 1.
}
call mw_sltermd (mw2, Memd[lterm+pdim2], Memd[lterm], pdim2)
ifnoerr (call mw_gwattrs (mw1, 1, "label", Memc[str], SZ_LINE))
call mw_swattrs (mw2, 1, "label", Memc[str])
ifnoerr (call mw_gwattrs (mw1, 1, "units", Memc[str], SZ_LINE))
call mw_swattrs (mw2, 1, "units", Memc[str])
ifnoerr (call mw_gwattrs (mw1, 1, "units_display", Memc[str], SZ_LINE))
call mw_swattrs (mw2, 1, "units_display", Memc[str])

# Write the MULTISPEC WCS
dtype = SMW_DTYPE(smw)
w1 = SMW_W1(smw)
dw = SMW_DW(smw)
nw = SMW_NW(smw)
z = SMW_Z(smw)
if (dtype == DCLOG) {
dw = log10 ((w1+(nw-1)*dw)/w1)/(nw-1)
w1 = log10 (w1)
}

call smwopn (mw2, smw, 0)
do i = 1, SMW_NSPEC(smw) {
ap = Memi[SMW_APS(smw)+i-1]
beam = Memi[SMW_BEAMS(smw)+i-1]
aplow = Memr[SMW_APLOW(smw)+2*i-2]
aphigh = Memr[SMW_APHIGH(smw)+2*i-2]
apid = Memi[SMW_APIDS(smw)+i-1]

call smw_mw (mw2, i, 1, mw, j, k)
call sprintf (Memc[key], SZ_FNAME, "spec%d")
call pargi (j)
call sprintf (Memc[str], SZ_LINE,
"%d %d %d %.14g %.14g %d %.14g %.2f %.2f")
call pargi (ap)
call pargi (beam)
call pargi (dtype)
call pargd (w1)
call pargd (dw)
call pargi (nw)
call pargd (z)
call pargr (aplow)
call pargr (aphigh)
call mw_swattrs (mw, 2, Memc[key], Memc[str])

if (apid != NULL)
call smw_sapid (mw2, i, 1, Memc[apid])

# This is used if there is a split MULTISPEC WCS.
if (SMW_APS(mw2) != NULL)
Memi[SMW_APS(mw2)+i-1] = ap
}

call smw_close (smw)
smw = mw2

call sfree (sp)
end
#@(#) smwgapid.x 1.1
include "smw.h"


# SMW_GAPID -- Get aperture id.

procedure smw_gapid (smw, index1, index2, apid, maxchar)

pointer smw #I SMW pointer
int index1 #I Spectrum index
int index2 #I Spectrum index
char apid[maxchar] #O Aperture id
int maxchar #I Maximum number of characters

pointer ptr

begin
switch (SMW_FORMAT(smw)) {
case SMW_ND:
call strcpy (Memc[SMW_APID(smw)], apid, maxchar)
case SMW_ES, SMW_MS:
if (index1 < 0 || index1 > SMW_NSPEC(smw))
call error (1, "smw_gapid: index out of bounds")

ptr = Memi[SMW_APIDS(smw)+index1-1]
if (index1 == 0 || ptr == NULL)
call strcpy (Memc[SMW_APID(smw)], apid, maxchar)
else
call strcpy (Memc[ptr], apid, maxchar)
}
end
#@(#) smwgwattrs.x 1.1
include
include "smw.h"


# SMW_GWATTRS -- Get spectrum attribute parameters.
# BE CAREFUL OF OUTPUT VARIABLES BEING THE SAME MEMORY ADDRESS!

procedure smw_gwattrs (smw, index1, index2, ap, beam, dtype, w1, dw, nw, z,
aplow, aphigh, coeff)

pointer smw # SMW pointer
int index1 # Spectrum index
int index2 # Spectrum index
int ap # Aperture number
int beam # Beam number
int dtype # Dispersion type
double w1 # Starting coordinate
double dw # Coordinate interval
int nw # Number of valid pixels
double z # Redshift factor
real aplow[2], aphigh[2] # Aperture limits
pointer coeff # Nonlinear coeff string (input/output)

int i, j, n, ip, sz_coeff, strlen(), ctoi(), ctor(), ctod()
double a, b
pointer sp, key, mw
errchk smw_mw, mw_gwattrs

data sz_coeff /SZ_LINE/

begin
call smark (sp)
call salloc (key, SZ_FNAME, TY_CHAR)

if (coeff == NULL)
call malloc (coeff, sz_coeff, TY_CHAR)
else
call realloc (coeff, sz_coeff, TY_CHAR)

# Determine parameters based on the SMW format.
switch (SMW_FORMAT(smw)) {
case SMW_ND:
call smw_mw (smw, index1, index2, mw, i, j)

dtype = SMW_DTYPE(smw)
nw = SMW_NW(smw)
w1 = SMW_W1(smw)
dw = SMW_DW(smw)
z = SMW_Z(smw)

ap = index1
beam = 0
aplow[1] = 1
aphigh[1] = 1
aplow[2] = 1
aphigh[2] = 1
if (SMW_LDIM(smw) > 1) {
aplow[1] = i - (SMW_NSUM(smw,1)-1) / 2
aphigh[1] = nint (aplow[1]) + SMW_NSUM(smw,1) - 1
aplow[1] = max (1, nint (aplow[1]))
aphigh[1] = min (SMW_LLEN(smw,2), nint (aphigh[1]))
}
if (SMW_LDIM(smw) > 2) {
aplow[2] = j - (SMW_NSUM(smw,2)-1) / 2
aphigh[2] = nint (aplow[2]) + SMW_NSUM(smw,2) - 1
aplow[2] = max (1, nint (aplow[2]))
aphigh[2] = min (SMW_LLEN(smw,3), nint (aphigh[2]))
}

Memc[coeff] = EOS
case SMW_ES:
call smw_mw (smw, index1, index2, mw, i, j)

dtype = SMW_DTYPE(smw)
nw = SMW_NW(smw)
w1 = SMW_W1(smw)
dw = SMW_DW(smw)
z = SMW_Z(smw)

ap = Memi[SMW_APS(smw)+index1-1]
beam = Memi[SMW_BEAMS(smw)+index1-1]
aplow[1] = Memr[SMW_APLOW(smw)+2*index1-2]
aphigh[1] = Memr[SMW_APHIGH(smw)+2*index1-2]
aplow[2] = Memr[SMW_APLOW(smw)+2*index1-1]
aphigh[2] = Memr[SMW_APHIGH(smw)+2*index1-1]

Memc[coeff] = EOS
case SMW_MS:
call smw_mw (smw, index1, index2, mw, i, j)

call sprintf (Memc[key], SZ_FNAME, "spec%d")
call pargi (i)

call mw_gwattrs (mw, 2, Memc[key], Memc[coeff], sz_coeff)
while (strlen (Memc[coeff]) == sz_coeff) {
sz_coeff = 2 * sz_coeff
call realloc (coeff, sz_coeff, TY_CHAR)
call mw_gwattrs (mw, 2, Memc[key], Memc[coeff], sz_coeff)
}

ip = 1
i = ctoi (Memc[coeff], ip, ap)
i = ctoi (Memc[coeff], ip, beam)
i = ctoi (Memc[coeff], ip, j)
i = ctod (Memc[coeff], ip, a)
i = ctod (Memc[coeff], ip, b)
i = ctoi (Memc[coeff], ip, n)
i = ctod (Memc[coeff], ip, z)
i = ctor (Memc[coeff], ip, aplow[1])
i = ctor (Memc[coeff], ip, aphigh[1])
aplow[2] = INDEF
aphigh[2] = INDEF
if (Memc[coeff+ip-1] != EOS)
call strcpy (Memc[coeff+ip], Memc[coeff], sz_coeff)
else
Memc[coeff] = EOS

if (j==DCLOG) {
if (abs(a)>20. || abs(a+(n-1)*b)>20.)
j = DCLINEAR
else {
a = 10D0 ** a
b = a * (10D0 ** ((n-1)*b) - 1) / (n - 1)
}
}

dtype = j
w1 = a
dw = b
nw = n
}

call sfree (sp)
end
#@(#) smwmerge.x 1.1
include
include "smw.h"


# SMW_MERGE -- Merge split MWCS array to a single MWCS.

procedure smw_merge (smw)

pointer smw #U Input split WCS, output single WCS

int i, pdim, naps, format, beam, dtype, dtype1, nw, nw1
int ap, axes[3]
double w1, dw, z, w11, dw1, z1
real aplow[2], aphigh[2]
pointer sp, key, val, term, coeff, mw, mw1, mw_open()
data axes/1,2,3/

begin
if (SMW_NMW(smw) == 1)
return

call smark (sp)
call salloc (key, SZ_FNAME, TY_CHAR)
call salloc (val, SZ_LINE, TY_CHAR)
call salloc (term, 15, TY_DOUBLE)
coeff = NULL

pdim = SMW_PDIM(smw)
naps = SMW_NSPEC(smw)
mw1 = SMW_MW(smw,0)

# Determine output WCS format.
format = SMW_ES
do i = 1, naps {
call smw_gwattrs (smw, i, 1, ap, beam, dtype, w1, dw, nw, z,
aplow, aphigh, coeff)
if (i == 1) {
dtype1 = dtype
w11 = w1
dw1 = dw
z1 = z
nw1 = nw
}
if (dtype>1||dtype!=dtype1||w1!=w11||dw!=dw1||nw!=nw1||z!=z1) {
format = SMW_MS
break
}
}

# Setup WCS.
switch (format) {
case SMW_ES:
mw = mw_open (NULL, pdim)
call mw_newsystem (mw, "equispec", pdim)
call mw_swtype (mw, axes, pdim, "linear", "")

case SMW_MS:
mw = mw_open (NULL, pdim)
call mw_newsystem (mw, "multispec", pdim)
call mw_swtype (mw, axes, pdim, "multispec", "")
if (pdim > 2)
call mw_swtype (mw, 3, 1, "linear", "")
}

ifnoerr (call mw_gwattrs (mw1, 1, "label", Memc[val], SZ_LINE))
call mw_swattrs (mw, 1, "label", Memc[val])
ifnoerr (call mw_gwattrs (mw1, 1, "units", Memc[val], SZ_LINE))
call mw_swattrs (mw, 1, "units", Memc[val])
ifnoerr (call mw_gwattrs (mw1, 1, "units_display", Memc[val], SZ_LINE))
call mw_swattrs (mw, 1, "units_display", Memc[val])
call mw_gltermd (mw1, Memd[term+pdim], Memd[term], pdim)
call mw_sltermd (mw, Memd[term+pdim], Memd[term], pdim)
call mw_gwtermd (mw1, Memd[term], Memd[term+pdim],
Memd[term+2*pdim], pdim)
Memd[term] = 1.
Memd[term+pdim] = w1 / (1 + z)
Memd[term+2*pdim] = dw / (1 + z)
call mw_swtermd (mw, Memd[term], Memd[term+pdim],
Memd[term+2*pdim], pdim)

# Set the SMW structure.
call smw_open (mw, smw, NULL)
do i = 1, naps {
call smw_gwattrs (smw, i, 1, ap, beam, dtype, w1, dw, nw, z,
aplow, aphigh, coeff)
call smw_swattrs (mw, i, 1, ap, beam, dtype, w1, dw, nw, z,
aplow, aphigh, coeff)
call smw_gapid (smw, i, 1, Memc[val], SZ_LINE)
call smw_sapid (mw, i, 1, Memc[val])
}

call smw_close (smw)
smw = mw

call mfree (coeff, TY_CHAR)
call sfree (sp)
end
#@(#) smwmultispec.x 1.1
include "smw.h"


# SMW_MULTISPEC -- Setup the MULTISPEC SMW parameters.

procedure smw_multispec (im, smw)

pointer im #I IMIO pointer
pointer smw #U MWCS pointer input SMW pointer output

int i, j, k
pointer sp, key, val, mw
errchk smw_open, smw_saxes, smw_sapid

begin
call smark (sp)
call salloc (key, SZ_FNAME, TY_CHAR)
call salloc (val, SZ_LINE, TY_CHAR)

call smw_open (smw, NULL, im)
do i = 1, SMW_NSPEC(smw) {
call smw_mw (smw, i, 1, mw, j, k)
call sprintf (Memc[key], SZ_FNAME, "APID%d")
call pargi (j)
ifnoerr (call imgstr (im, Memc[key], Memc[val], SZ_LINE))
call smw_sapid (smw, i, 1, Memc[val])
}

call sfree (sp)
end
#@(#) smwmw.x 1.1
include "smw.h"


# SMW_MW -- Get MWCS pointer and coordinates from spectrum line and band

procedure smw_mw (smw, line, band, mw, x, y)

pointer smw #I SMW pointer
int line #I Spectrum line
int band #I Spectrum band
pointer mw #O MWCS pointer
int x, y #O MWCS coordinates

real mw_c1tranr()

begin
if (line < 1 || line > SMW_NSPEC(smw))
call error (1, "smw_mw: spectrum not found")

switch (SMW_FORMAT(smw)) {
case SMW_ND:
mw = SMW_MW(smw,0)
x = mod (line - 1, SMW_LLEN(smw,2)) + 1
y = (line - 1) / SMW_LLEN(smw,2) + band
default:
if (SMW_NMW(smw) == 1) {
mw = SMW_MW(smw,0)
x = line
y = band
if (SMW_CTLP(smw) != NULL)
x = nint (mw_c1tranr (SMW_CTLP(smw), real(line)))
} else {
mw = SMW_MW(smw,(line-1)/SMW_NSPLIT)
x = mod (line - 1, SMW_NSPLIT) + 1
y = band
}
}
end
#@(#) smwnd.x 1.1
include
include "smw.h"


# SMW_ND -- Setup the NDSPEC SMW.
# If there is only one spectrum convert it to EQUISPEC if possible.

procedure smw_nd (im, smw)

pointer im #I IMIO pointer
pointer smw #U MWCS pointer input SMW pointer output

errchk smw_open, smw_daxis, smw_ndes

begin
call smw_open (smw, NULL, im)
if (SMW_NSPEC(smw) == 1)
call smw_ndes (im, smw)
end
#@(#) smwndes.x 1.1
include
include "smw.h"


# SMW_NDES -- Convert NDSPEC WCS into EQUISPEC WCS.
# This requires that the logical dispersion axis be 1.

procedure smw_ndes (im, smw)

pointer im #I IMIO pointer
pointer smw #U Input NDSPEC SMW, output EQUISPEC SMW

int i, pdim1, pdim2, daxis, ap, beam, dtype, nw, axes[2]
real aplow[2], aphigh[2]
double w1, dw, z
pointer sp, key, str, lterm1, lterm2, coeff, mw1, mw2, mw_open()
errchk mw_open, mw_gltermd, mw_gwtermd, smw_open, smw_saxes, smw_gwattrs
data axes/1,2/, coeff/NULL/

begin
# Require the dispersion to be along the first logical axis.
if (SMW_LAXIS(smw,1) != 1)
return

call smark (sp)
call salloc (key, SZ_FNAME, TY_CHAR)
call salloc (str, SZ_LINE, TY_CHAR)
call salloc (lterm1, 15, TY_DOUBLE)
call salloc (lterm2, 15, TY_DOUBLE)

# Set the MWCS. Only the logical and world transformations along
# the dispersion axis are transfered.

pdim1 = SMW_PDIM(smw)
pdim2 = IM_NDIM(im)
daxis = SMW_PAXIS(smw,1)
mw1 = SMW_MW(smw,0)

mw2 = mw_open (NULL, pdim2)
call mw_newsystem (mw2, "equispec", pdim2)
call mw_swtype (mw2, axes, pdim2, "linear", "")
ifnoerr (call mw_gwattrs (mw1, daxis, "label", Memc[str], SZ_LINE))
call mw_swattrs (mw2, 1, "label", Memc[str])
ifnoerr (call mw_gwattrs (mw1, daxis, "units", Memc[str], SZ_LINE))
call mw_swattrs (mw2, 1, "units", Memc[str])
ifnoerr (call mw_gwattrs (mw1, daxis, "units_display", Memc[str],
SZ_LINE))
call mw_swattrs (mw2, 1, "units_display", Memc[str])

call mw_gltermd (mw1, Memd[lterm1+pdim1], Memd[lterm1], pdim1)
call mw_gltermd (mw2, Memd[lterm2+pdim2], Memd[lterm2], pdim2)
Memd[lterm2] = Memd[lterm1+daxis-1]
Memd[lterm2+pdim2] = Memd[lterm1+pdim1+(pdim1+1)*(daxis-1)]
call mw_sltermd (mw2, Memd[lterm2+pdim2], Memd[lterm2], pdim2)

call mw_gwtermd (mw1, Memd[lterm1], Memd[lterm1+pdim1],
Memd[lterm1+2*pdim1], pdim1)
call mw_gwtermd (mw2, Memd[lterm2], Memd[lterm2+pdim2],
Memd[lterm2+2*pdim2], pdim2)
Memd[lterm2] = Memd[lterm1+daxis-1]
Memd[lterm2+pdim2] = Memd[lterm1+pdim1+daxis-1]
Memd[lterm2+2*pdim2] = Memd[lterm1+2*pdim1+(pdim1+1)*(daxis-1)]
call mw_swtermd (mw2, Memd[lterm2], Memd[lterm2+pdim2],
Memd[lterm2+2*pdim2], pdim2)

# Set the EQUISPEC SMW.
IM_LEN(im,2) = SMW_NSPEC(smw)
IM_LEN(im,3) = SMW_NBANDS(smw)
call smw_open (mw2, NULL, im)
do i = 1, SMW_NSPEC(smw) {
call smw_gwattrs (smw, i, 1, ap, beam, dtype, w1, dw, nw, z,
aplow, aphigh, coeff)
call smw_swattrs (mw2, i, 1, ap, beam, dtype, w1, dw, nw, z,
aplow, aphigh, Memc[coeff])
}
call mfree (coeff, TY_CHAR)

call smw_close (smw)
smw = mw2

call sfree (sp)
end
#@(#) smwnewcopy.x 1.1
include "smw.h"


# SMW_NEWCOPY -- Make a new copy of an SMW structure.

pointer procedure smw_newcopy (smw)

pointer smw #I SMW pointer to copy
pointer smwcopy #O SMW copy

int i
pointer mw_newcopy()

begin
smwcopy = SMW_MW(smw,0)
call smw_open (smwcopy, smw, NULL)

if (SMW_APS(smwcopy) != NULL)
call amovi (Memi[SMW_APS(smw)], Memi[SMW_APS(smwcopy)], i)

do i = 0, SMW_NMW(smw)-1
SMW_MW(smwcopy,i) = mw_newcopy (SMW_MW(smw,i))

return (smwcopy)
end
#@(#) smwoldms.x 1.1
include
include "smw.h"


# SMW_OLDMS -- Convert old multispec format into MULTISPEC SMW.

procedure smw_oldms (im, smw)

pointer im #I IMIO pointer
pointer smw #U Input MWCS pointer, output SMW pointer

int i, j, k, nchar, ap, beam, dtype, nw, axes[2]
double w1, dw, z
real aplow[2], aphigh[2]
pointer sp, key, val, lterm, mw, mw_open()
int imgeti(), mw_stati(), ctoi(), ctor(), ctod(), imofnlu(), imgnfn()
errchk imgstr, mw_gltermd, mw_sltermd
data axes/1,2/

begin
call smark (sp)
call salloc (key, SZ_FNAME, TY_CHAR)
call salloc (val, SZ_LINE, TY_CHAR)
call salloc (lterm, 12, TY_DOUBLE)

# Set the basic multispec MWCS
i = mw_stati (smw, MW_NDIM)
j = max (2, i)
mw = mw_open (NULL, j)
call mw_newsystem (mw, "multispec", j)
call mw_swtype (mw, axes, 2, "multispec", "")
if (j > 2)
call mw_swtype (mw, 3, 1, "linear", "")
call mw_gltermd (smw, Memd[lterm+j], Memd[lterm], i)
if (i == 1) {
Memd[lterm+1] = 0.
Memd[lterm+3] = 0.
Memd[lterm+4] = 0.
Memd[lterm+5] = 1.
}
call mw_sltermd (mw, Memd[lterm+j], Memd[lterm], j)

iferr (dtype = imgeti (im, "DC-FLAG"))
dtype = -1
else {
call mw_swattrs (mw, 1, "label", "Wavelength")
call mw_swattrs (mw, 1, "units", "Angstroms")
}

call mw_close (smw)
smw = mw

# Set the SMW data structure.
call smw_open (smw, NULL, im)
do i = 1, SMW_NSPEC(smw) {
call smw_mw (smw, i, 1, mw, j, k)
call sprintf (Memc[key], SZ_FNAME, "APNUM%d")
call pargi (j)
call imgstr (im, Memc[key], Memc[val], SZ_LINE)
call imdelf (im, Memc[key])

k = 1
nchar = ctoi (Memc[val], k, ap)
nchar = ctoi (Memc[val], k, beam)
nchar = ctod (Memc[val], k, w1)
nchar = ctod (Memc[val], k, dw)
nchar = ctoi (Memc[val], k, nw)
if (ctor (Memc[val], k, aplow[1]) == 0)
aplow[1] = INDEF
if (ctor (Memc[val], k, aphigh[1]) == 0)
aphigh[1] = INDEF
z = 0.

k = dtype
if (k==1 && (abs(w1)>20. || abs(w1+(nw-1)*dw)>20.))
k = 0
call smw_swattrs (smw, i, 1, ap, beam, k, w1, dw, nw, z,
aplow, aphigh, "")

call sprintf (Memc[key], SZ_FNAME, "APID%d")
call pargi (j)
ifnoerr (call imgstr (im, Memc[key], Memc[val], SZ_LINE)) {
call smw_sapid (smw, i, 1, Memc[val])
call imdelf (im, Memc[key])
}
}

# Delete old parameters
i = imofnlu (im,
"DISPAXIS,APFORMAT,BEAM-NUM,DC-FLAG,W0,WPC,NP1,NP2")
while (imgnfn (i, Memc[key], SZ_FNAME) != EOF)
call imdelf (im, Memc[key])
call imcfnl (i)

# Update MWCS
call smw_saveim (smw, im)

call sfree (sp)
end
#@(#) smwonedspec.x 1.1
include
include "smw.h"


# SMW_ONEDSPEC -- Convert old "onedspec" format to EQUISPEC.

procedure smw_onedspec (im, smw)

pointer im #I IMIO pointer
pointer smw #U MWCS pointer input SMW pointer output

int i, dtype, ap, beam, nw, imgeti(), imofnlu(), imgnfn()
real aplow[2], aphigh[2], imgetr(), mw_c1tranr()
double ltm, ltv, r, w, dw, z, imgetd()
pointer sp, key, mw, ct, mw_openim(), mw_sctran()
errchk smw_open, smw_saxes, mw_gwtermd, mw_sctran

begin
call smark (sp)
call salloc (key, SZ_FNAME, TY_CHAR)

# Convert old W0/WPC keywords if needed.
mw = smw
iferr (w = imgetd (im, "CRVAL1")) {
ifnoerr (w = imgetd (im, "W0")) {
dw = imgetd (im, "WPC")
iferr (ltm = imgetd (im, "LTM1_1"))
ltm = 1
iferr (ltv = imgetd (im, "LTV1"))
ltv = 0
r = ltm + ltv
dw = dw / ltm
call imaddd (im, "CRPIX1", r)
call imaddd (im, "CRVAL1", w)
call imaddd (im, "CD1_1", dw)
call imaddd (im, "CDELT1", dw)
call mw_close(mw)
mw = mw_openim (im)
}
}

# Get dispersion and determine number of valid pixels.
call mw_gwtermd (mw, r, w, dw, 1)
w = w - (r - 1) * dw
r = 1
call mw_swtermd (mw, r, w, dw, 1)
ct = mw_sctran (mw, "logical", "physical", 1)
nw = max (mw_c1tranr (ct, 1.), mw_c1tranr (ct, real (IM_LEN(im,1))))
call mw_ctfree (ct)

iferr (dtype = imgeti (im, "DC-FLAG"))
dtype = DCNO
if (dtype==DCLOG) {
if (abs(w)>20. || abs(w+(nw-1)*dw)>20.)
dtype = DCLINEAR
else {
w = 10D0 ** w
dw = w * (10D0 ** ((nw-1)*dw) - 1) / (nw - 1)
}
}

# Convert to EQUISPEC system.
call mw_swattrs (mw, 0, "system", "equispec")
if (dtype != -1) {
call mw_swattrs (mw, 1, "label", "Wavelength")
call mw_swattrs (mw, 1, "units", "angstroms")
}

# Set the SMW data structure.
call smw_open (smw, NULL, im)

# Determine the aperture parameters.
iferr (beam = imgeti (im, "BEAM-NUM"))
beam = 1
iferr (ap = imgeti (im, "APNUM"))
ap = beam
iferr (aplow[1] = imgetr (im, "APLOW"))
aplow[1] = INDEF
iferr (aphigh[1] = imgetr (im, "APHIGH"))
aphigh[1] = INDEF
iferr (z = imgetd (im, "DOPCOR"))
z = 0.

call smw_swattrs (smw, 1, 1, ap, beam, dtype, w, dw, nw, z,
aplow, aphigh, "")

# Delete old parameters
i = imofnlu (im,
"BEAM-NUM,APNUM,APLOW,APHIGH,DOPCOR,DC-FLAG,W0,WPC,NP1,NP2")
while (imgnfn (i, Memc[key], SZ_FNAME) != EOF)
call imdelf (im, Memc[key])
call imcfnl (i)

# Update MWCS
call smw_saveim (smw, im)

call sfree (sp)
end
#@(#) smwopen.x 1.1
include "smw.h"


# SMW_OPEN -- Open SMW structure.
# The basic MWCS pointer and a template SMW pointer or image is input
# and the SMW pointer is returned in its place.

procedure smw_open (mw, smw1, im)

pointer mw #U MWCS pointer input and SMW pointer output
pointer smw1 #I Template SMW pointer
pointer im #I Template IMIO pointer

int i, nspec, nmw, format, strdic()
pointer sp, sys, smw, mw_sctran(), mw_newcopy()
errchk smw_daxis, smw_saxes, mw_sctran

begin
call smark (sp)
call salloc (sys, SZ_FNAME, TY_CHAR)

call mw_gwattrs (mw, 0, "system", Memc[sys], SZ_FNAME)
format = strdic (Memc[sys], Memc[sys], SZ_FNAME, SMW_FORMATS)

call calloc (smw, SMW_LEN(1), TY_STRUCT)
call malloc (SMW_APID(smw), SZ_LINE, TY_CHAR)
SMW_FORMAT(smw) = format
SMW_DTYPE(smw) = INDEFI
SMW_NMW(smw) = 1
SMW_MW(smw,0) = mw

switch (format) {
case SMW_ND:
call smw_daxis (smw, im, INDEFI, INDEFI, INDEFI)
call smw_saxes (smw, smw1, im)

case SMW_ES:
call smw_saxes (smw, smw1, im)

nspec = SMW_NSPEC(smw)
call malloc (SMW_APS(smw), nspec, TY_INT)
call malloc (SMW_BEAMS(smw), nspec, TY_INT)
call malloc (SMW_APLOW(smw), 2*nspec, TY_REAL)
call malloc (SMW_APHIGH(smw), 2*nspec, TY_REAL)
call calloc (SMW_APIDS(smw), nspec, TY_POINTER)
if (SMW_PDIM(smw) > 1)
SMW_CTLP(smw) = mw_sctran (mw, "logical", "physical", 2)

case SMW_MS:
call smw_saxes (smw, smw1, im)

nspec = SMW_NSPEC(smw)
call calloc (SMW_APIDS(smw), nspec, TY_POINTER)
if (SMW_PDIM(smw) > 1)
SMW_CTLP(smw) = mw_sctran (mw, "logical", "physical", 2)

nmw = 1 + (nspec - 1) / SMW_NSPLIT
if (nmw > 1) {
call realloc (smw, SMW_LEN(nmw), TY_STRUCT)
call calloc (SMW_APS(smw), nspec, TY_INT)
}
do i = 1, nmw-1
SMW_MW(smw,i) = mw_newcopy (mw)
SMW_NMW(smw) = nmw
}

mw = smw

call sfree (sp)
end
#@(#) smwopenim.x 1.1
include
include
include

define SYSTEMS "|equispec|multispec|physical|image|world|linear|"


# SMW_OPENIM -- Open the spectral MWCS for various input formats.

pointer procedure smw_openim (im)

pointer im #I Image pointer
pointer mw #O MWCS pointer

pointer sp, system, mw_openim()
bool streq()
int i, wcsdim, sys, strdic(), mw_stati()
errchk mw_openim, smw_oldms, smw_linear

begin
call smark (sp)
call salloc (system, SZ_FNAME, TY_CHAR)

# Workaround for truncation of header during image header copy.
IM_HDRLEN(im) = IM_LENHDRMEM(im)

# Force higher dimensions to length 1.
do i = IM_NDIM(im) + 1, 3
IM_LEN(im,i) = 1

mw = mw_openim (im)
call mw_seti (mw, MW_USEAXMAP, NO)
wcsdim = mw_stati (mw, MW_NDIM)
call mw_gwattrs (mw, 0, "system", Memc[system], SZ_FNAME)
sys = strdic (Memc[system], Memc[system], SZ_FNAME, SYSTEMS)

# Set various input systems.
switch (sys) {
case 1:
call smw_equispec (im, mw)
case 2:
call smw_multispec (im, mw)
default:
if (sys == 0) {
call eprintf (
"WARNING: Unknown coordinate system `%s' - assuming `linear'.\n")
call pargstr (Memc[system])
} else if (sys == 3)
call mw_newsystem (mw, "image", wcsdim)

# Old "multispec" format.
ifnoerr (call imgstr (im, "APFORMAT", Memc[system], SZ_FNAME)) {
if (streq (Memc[system], "onedspec"))
call smw_onedspec (im, mw)
else
call smw_oldms (im, mw)

# Old "onedspec" format or other 1D image.
} else if (wcsdim == 1) {
call smw_onedspec (im, mw)

# N-dimensional image.
} else
call smw_nd (im, mw)
}

call sfree (sp)
return (mw)
end
#@(#) smwsapid.x 1.1
include "smw.h"


# SMW_SAPID -- Set aperture id.

procedure smw_sapid (smw, index1, index2, apid)

pointer smw #I SMW pointer
int index1 #I Spectrum index
int index2 #I Spectrum index
char apid[ARB] #I Aperture id

pointer ptr
bool streq()
errchk malloc

begin
switch (SMW_FORMAT(smw)) {
case SMW_ND:
call strcpy (apid, Memc[SMW_APID(smw)], SZ_LINE)
case SMW_ES, SMW_MS:
if (index1 < 0 || index1 > SMW_NSPEC(smw))
call error (1, "smw_sapid: index out of bounds")

if (index1 == 0)
call strcpy (apid, Memc[SMW_APID(smw)], SZ_LINE)

else {
ptr = Memi[SMW_APIDS(smw)+index1-1]
if (streq (apid, Memc[SMW_APID(smw)]))
call mfree (ptr, TY_CHAR)
else {
if (ptr == NULL)
call malloc (ptr, SZ_LINE, TY_CHAR)
call strcpy (apid, Memc[ptr], SZ_LINE)
}
Memi[SMW_APIDS(smw)+index1-1] = ptr
}
}
end
#@(#) smwsaveim.x 1.1
include
include
include "smw.h"


# SMW_SAVEIM -- Save spectral WCS in image header.
# The input and output formats are EQUISPEC and MULTISPEC. A split input
# MULTISPEC WCS is first merged to a single EQUISPEC or MULTISPEC WCS.
# An input MULTISPEC WCS is converted to EQUISPEC output if possible.

procedure smw_saveim (smw, im)

pointer smw # SMW pointer
pointer im # Image pointer

int i, j, format, nl, pdim, pdim1, beam, dtype, dtype1, nw, nw1
int ap, axes[3]
real aplow[2], aphigh[2]
double v, m, w1, dw, z, w11, dw1, z1
pointer sp, key, str1, str2, axmap, lterm, coeff, mw, mw1

bool strne(), fp_equald()
int imaccf()
pointer mw_open()
errchk smw_merge, smw_c2trand, imdelf
data axes/1,2,3/

begin
call smark (sp)
call salloc (key, SZ_FNAME, TY_CHAR)
call salloc (str1, SZ_LINE, TY_CHAR)
call salloc (str2, SZ_LINE, TY_CHAR)
call salloc (axmap, 6, TY_INT)
call salloc (lterm, 15, TY_DOUBLE)
coeff = NULL

# Merge split WCS into a single WCS.
call smw_merge (smw)

mw = SMW_MW(smw,0)
pdim = SMW_PDIM(smw)
format = SMW_FORMAT(smw)
if (IM_NDIM(im) == 1)
nl = 1
else
nl = IM_LEN(im,2)

# Check if MULTISPEC WCS can be converted to EQUISPEC.
if (format == SMW_MS) {
format = SMW_ES
do i = 1, nl {
call smw_gwattrs (smw, i, 1, ap, beam, dtype, w1, dw, nw, z,
aplow, aphigh, coeff)
if (i == 1) {
dtype1 = dtype
w11 = w1
dw1 = dw
z1 = z
nw1 = nw
}
if (dtype>1||dtype!=dtype1||!fp_equald(w1,w11)||
!fp_equald(dw,dw1)||nw!=nw1||!fp_equald(z,z1)) {
format = SMW_MS
break
}
}
}

# Save WCS in desired format.
switch (format) {
case SMW_ND:
if (SMW_DTYPE(smw) != -1)
call imaddi (im, "DC-FLAG", SMW_DTYPE(smw))
else if (imaccf (im, "DC-FLAG") == YES)
call imdelf (im, "DC-FLAG")
if (imaccf (im, "DISPAXIS") == YES)
call imaddi (im, "DISPAXIS", SMW_PAXIS(smw,1))

call smw_gapid (smw, 1, 1, IM_TITLE(im), SZ_IMTITLE)
call mw_saveim (mw, im)

case SMW_ES:
# Save aperture information.
do i = 1, nl {
call smw_gwattrs (smw, i, 1, ap, beam, dtype, w1, dw, nw, z,
aplow, aphigh, coeff)
call sprintf (Memc[key], SZ_FNAME, "APNUM%d")
call pargi (i)
call sprintf (Memc[str1], SZ_LINE, "%d %d")
call pargi (ap)
call pargi (beam)
if (!IS_INDEF(aplow[1]) || !IS_INDEF(aphigh[1])) {
call sprintf (Memc[str2], SZ_LINE, " %.2f %.2f")
call pargr (aplow[1])
call pargr (aphigh[1])
call strcat (Memc[str2], Memc[str1], SZ_LINE)
if (!IS_INDEF(aplow[2]) || !IS_INDEF(aphigh[2])) {
call sprintf (Memc[str2], SZ_LINE, " %.2f %.2f")
call pargr (aplow[2])
call pargr (aphigh[2])
call strcat (Memc[str2], Memc[str1], SZ_LINE)
}
}
call imastr (im, Memc[key], Memc[str1])
if (i == 1) {
iferr (call imdelf (im, "APID1"))
;
}
call smw_gapid (smw, i, 1, Memc[str1], SZ_LINE)
if (Memc[str1] != EOS) {
if (strne (Memc[str1], IM_TITLE(im))) {
if (nl == 1) {
call imastr (im, "MSTITLE", IM_TITLE(im))
call strcpy (Memc[str1], IM_TITLE(im), SZ_IMTITLE)
} else {
call sprintf (Memc[key], SZ_FNAME, "APID%d")
call pargi (i)
call imastr (im, Memc[key], Memc[str1])
}
}
}
}

# Delete unnecessary aperture information.
do i = nl+1, ARB {
call sprintf (Memc[key], SZ_FNAME, "APNUM%d")
call pargi (i)
iferr (call imdelf (im, Memc[key]))
break
call sprintf (Memc[key], SZ_FNAME, "APID%d")
call pargi (i)
iferr (call imdelf (im, Memc[key]))
;
}

# Add dispersion parameters to image.
if (dtype != -1)
call imaddi (im, "DC-FLAG", dtype)
else if (imaccf (im, "DC-FLAG") == YES)
call imdelf (im, "DC-FLAG")
if (nw < IM_LEN(im,1))
call imaddi (im, "NP2", nw)
else if (imaccf (im, "NP2") == YES)
call imdelf (im, "NP2")

# Setup EQUISPEC WCS.
pdim1 = max (IM_NDIM(im), IM_NPHYSDIM(im))
mw1 = mw_open (NULL, pdim1)
call mw_newsystem (mw1, "equispec", pdim1)
call mw_swtype (mw1, axes, pdim1, "linear", "")
ifnoerr (call mw_gwattrs (mw, 1, "label", Memc[str1], SZ_LINE))
call mw_swattrs (mw1, 1, "label", Memc[str1])
ifnoerr (call mw_gwattrs (mw, 1, "units", Memc[str1], SZ_LINE))
call mw_swattrs (mw1, 1, "units", Memc[str1])
ifnoerr (call mw_gwattrs (mw, 1, "units_display", Memc[str1],
SZ_LINE))
call mw_swattrs (mw1, 1, "units_display", Memc[str1])
call mw_gltermd (mw, Memd[lterm+pdim], Memd[lterm], pdim)
v = Memd[lterm]
m = Memd[lterm+pdim]
call mw_gltermd (mw1, Memd[lterm+pdim1], Memd[lterm], pdim1)
Memd[lterm] = v
Memd[lterm+pdim1] = m
call mw_sltermd (mw1, Memd[lterm+pdim1], Memd[lterm], pdim1)
call mw_gwtermd (mw1, Memd[lterm], Memd[lterm+pdim1],
Memd[lterm+2*pdim1], pdim1)
Memd[lterm] = 1.
w1 = w1 / (1 + z)
dw = dw / (1 + z)
if (dtype == DCLOG) {
dw = log10 ((w1 + (nw - 1) * dw) / w1) / (nw - 1)
w1 = log10 (w1)
}
Memd[lterm+pdim1] = w1
Memd[lterm+2*pdim1] = dw
call mw_swtermd (mw1, Memd[lterm], Memd[lterm+pdim1],
Memd[lterm+2*pdim1], pdim1)
call mw_saveim (mw1, im)
call mw_close (mw1)

case SMW_MS:
# Delete any APNUM keywords. If there is only one spectrum
# define the axis mapping.

do j = 1, ARB {
call sprintf (Memc[key], SZ_FNAME, "APNUM%d")
call pargi (j)
iferr (call imdelf (im, Memc[key]))
break
}
if (IM_NDIM(im) == 1) {
call aclri (Memi[axmap], 2*pdim)
Memi[axmap] = 1
call mw_saxmap (mw, Memi[axmap], Memi[axmap+pdim], pdim)
}

# Set aperture ids.
do i = 1, nl {
if (i == 1) {
iferr (call imdelf (im, "APID1"))
;
}
call smw_gapid (smw, i, 1, Memc[str1], SZ_LINE)
if (Memc[str1] != EOS) {
if (strne (Memc[str1], IM_TITLE(im))) {
if (nl == 1) {
call imastr (im, "MSTITLE", IM_TITLE(im))
call strcpy (Memc[str1], IM_TITLE(im), SZ_IMTITLE)
} else {
call sprintf (Memc[key], SZ_FNAME, "APID%d")
call pargi (i)
call imastr (im, Memc[key], Memc[str1])
}
}
}
}

do i = nl+1, ARB {
call sprintf (Memc[key], SZ_FNAME, "APID%d")
call pargi (i)
iferr (call imdelf (im, Memc[key]))
break
}

call mw_saveim (mw, im)
}

call mfree (coeff, TY_CHAR)
call sfree (sp)
end
#@(#) smwsaxes.x 1.1
include
include
include "smw.h"


# SMW_SAXES -- Set axes parameters based on previously set dispersion axis.
# If the dispersion axis has been excluded for NDSPEC allow another axis to
# be chosen with a warning. For EQUISPEC and MULTISPEC require the dispersion
# to be 1 and also to be present.

procedure smw_saxes (smw, smw1, im)

pointer smw #I SMW pointer
pointer smw1 #I Template SMW pointer
pointer im #I Template IMIO pointer

int i, pdim, ldim, paxis, laxis, nw, dtype, nspec
real smw_c1tranr()
double w1, dw
pointer sp, axno, axval, r, w, cd, mw, ct, smw_sctran()
int mw_stati(), imgeti()

begin
# If a template SMW pointer is specified just copy the axes parameters.
if (smw1 != NULL) {
call strcpy (Memc[SMW_APID(smw1)], Memc[SMW_APID(smw)], SZ_LINE)
SMW_NSPEC(smw) = SMW_NSPEC(smw1)
SMW_NBANDS(smw) = SMW_NBANDS(smw1)
SMW_TRANS(smw) = SMW_TRANS(smw1)
call amovi (SMW_PAXIS(smw1,1), SMW_PAXIS(smw,1), 3)
SMW_LDIM(smw) = SMW_LDIM(smw1)
call amovi (SMW_LAXIS(smw1,1), SMW_LAXIS(smw,1), 3)
call amovi (SMW_LLEN(smw1,1), SMW_LLEN(smw,1), 3)
call amovi (SMW_NSUM(smw1,1), SMW_NSUM(smw,1), 2)

mw = SMW_MW(smw,0)
SMW_PDIM(smw) = mw_stati (mw, MW_NDIM)
if (SMW_PDIM(smw) > SMW_PDIM(smw1))
do i = SMW_PDIM(smw1)+1, SMW_PDIM(smw)
SMW_PAXIS(smw,i) = i

return
}

call smark (sp)
call salloc (axno, 3, TY_INT)
call salloc (axval, 3, TY_INT)
call aclri (Memi[axno], 3)

# Determine the dimensions.
mw = SMW_MW(smw,0)
pdim = mw_stati (mw, MW_NDIM)
ldim = IM_NDIM(im)
call mw_gaxmap (mw, Memi[axno], Memi[axval], pdim)

# Set the physical dispersion axis.
switch (SMW_FORMAT(smw)) {
case SMW_ND:
call salloc (r, pdim, TY_DOUBLE)
call salloc (w, pdim, TY_DOUBLE)
call salloc (cd, pdim*pdim, TY_DOUBLE)

# Check for a transposed 2D image.
SMW_TRANS(smw) = NO
if (pdim == 2) {
call mw_gltermd (mw, Memd[cd], Memd[w], pdim)
if (Memd[cd] == 0D0 && Memd[cd+3] == 0D0) {
Memd[cd] = Memd[cd+2]
Memd[cd+2] = 0.
Memd[cd+3] = Memd[cd+1]
Memd[cd+1] = 0.
call mw_sltermd (mw, Memd[cd], Memd[w], pdim)
paxis = SMW_PAXIS(smw,1)
if (paxis == 1)
SMW_PAXIS(smw,1) = 2
else
SMW_PAXIS(smw,1) = 1
SMW_TRANS(smw) = YES
}
}

# If the dispersion axis is of length 1 or has been excluded find
# the first longer axis and print a warning.

paxis = SMW_PAXIS(smw,1)
i = max (1, min (pdim, paxis))
laxis = max (1, Memi[axno+i-1])
if (IM_LEN(im,laxis) == 1)
do laxis = 1, ldim
if (IM_LEN(im,laxis) != 1)
break

# Determine the number of spectra.
nspec = 1
do i = 1, ldim
if (i != laxis)
nspec = nspec * IM_LEN(im,i)
SMW_NSPEC(smw) = nspec
SMW_NBANDS(smw) = 1

i = paxis
do paxis = 1, pdim
if (Memi[axno+paxis-1] == laxis)
break

if (i != paxis && nspec > 1) {
call eprintf (
"WARNING: Dispersion axis %d not found. Using axis %d.\n")
call pargi (i)
call pargi (paxis)
}

# Set the dispersion system.
call mw_gwtermd (mw, Memd[r], Memd[w], Memd[cd], pdim)
if (SMW_TRANS(smw) == YES) {
w1 = Memd[r]
Memd[r] = Memd[r+1]
Memd[r+1] = w1
Memd[cd] = Memd[cd+1]
Memd[cd+1] = 0.
Memd[cd+3] = Memd[cd+2]
Memd[cd+2] = 0.
}
do i = 0, pdim-1 {
dw = Memd[cd+i*(pdim+1)]
if (dw == 0D0)
Memd[cd+i*(pdim+1)] = 1D0
}
call mw_swtermd (mw, Memd[r], Memd[w], Memd[cd], pdim)

dw = Memd[cd+(paxis-1)*(pdim+1)]
w1 = Memd[w+paxis-1] - (Memd[r+paxis-1] - 1) * dw
nw = IM_LEN(im,laxis)

i = 2 ** (paxis - 1)
ct = smw_sctran (smw, "logical", "physical", i)
nw = max (smw_c1tranr (ct, 0.5), smw_c1tranr (ct, nw+0.5))
call smw_ctfree (ct)

iferr (dtype = imgeti (im, "DC-FLAG"))
dtype = DCNO
if (dtype==DCLOG) {
if (abs(w1)>20. || abs(w1+(nw-1)*dw)>20.)
dtype = DCLINEAR
else {
w1 = 10D0 ** w1
dw = w1 * (10D0 ** ((nw-1)*dw) - 1) / (nw - 1)
}
}

SMW_DTYPE(smw) = INDEFI
call smw_swattrs (smw, 1, 1, INDEFI, INDEFI,
dtype, w1, dw, nw, 0D0, INDEFR, INDEFR, "")
case SMW_ES, SMW_MS:
paxis = 1
i = Memi[axno+1]
if (i == 0)
SMW_NSPEC(smw) = 1
else
SMW_NSPEC(smw) = IM_LEN(im,i)
i = Memi[axno+2]
if (i == 0)
SMW_NBANDS(smw) = 1
else
SMW_NBANDS(smw) = IM_LEN(im,i)
}

# Check and set the physical and logical dispersion axes.
laxis = Memi[axno+paxis-1]
if (laxis == 0) {
if (Memi[axval+paxis-1] == 0)
laxis = paxis
else
call error (1, "No dispersion axis")
}

SMW_PDIM(smw) = pdim
SMW_LDIM(smw) = ldim
SMW_PAXIS(smw,1) = paxis
SMW_LAXIS(smw,1) = laxis
SMW_LLEN(smw,1) = IM_LEN(im,laxis)
SMW_LLEN(smw,2) = 1
SMW_LLEN(smw,3) = 1

# Set the spatial axes.
i = 2
do laxis = 1, ldim {
if (laxis != SMW_LAXIS(smw,1)) {
do paxis = 1, pdim
if (Memi[axno+paxis-1] == laxis)
break
SMW_PAXIS(smw,i) = paxis
SMW_LAXIS(smw,i) = laxis
SMW_LLEN(smw,i) = IM_LEN(im,laxis)
i = i + 1
}
}

# Set the default title.
call smw_sapid (smw, 0, 1, IM_TITLE(im))

call sfree (sp)
end
#@(#) smwsctran.x 1.1
include "smw.h"


# SMW_SCTRAN -- Set the SMW coordinate system transformation.
# This routine sets up the SMW_CT structure which may consist of multiple
# pointers if the MWCS is split. If the MWCS is not split then the MWCS
# transformation routine is used directly. However if the MWCS is split then
# there may need to be an intermediate mapping from the input coordinate to
# the physical coordinate in which the split MWCS is defined as well as a
# transformation for each WCS piece.

pointer procedure smw_sctran (smw, system1, system2, axbits)

pointer smw #I SMW pointer
char system1[ARB] #I Input coordinate system
char system2[ARB] #I Output coordinate system
int axbits #I Bitmap defining axes to be transformed
pointer ct #O SMW CT pointer

int i, cttype, nct, strdic()
pointer mw_sctran()
errchk mw_sctran

begin
# Determine the coordinate transformation type and setup the structure.
cttype = 10 * (strdic (system1, system1, ARB, SMW_CTTYPES) + 1) +
strdic (system2, system2, ARB, SMW_CTTYPES) + 1

nct = SMW_NMW(smw)
if (cttype == SMW_LP || cttype == SMW_PL)
nct = 1

call calloc (ct, SMW_CTLEN(nct), TY_STRUCT)
SMW_SMW(ct) = smw
SMW_CTTYPE(ct) = cttype
SMW_NCT(ct) = nct

# If the MWCS is not split use the MWCS transformation directly.
if (nct == 1) {
SMW_CT(ct,0) = mw_sctran (SMW_MW(smw,0), system1, system2, axbits)
return(ct)
}

# If there is a split MWCS then setup the intermediary transformation.
switch (cttype) {
case SMW_LW:
SMW_CTL(ct) = mw_sctran (SMW_MW(smw,0), system1, "physical",
axbits)
do i = 0, nct-1
SMW_CT(ct,i) = mw_sctran (SMW_MW(smw,i), "physical",
system2, axbits)
case SMW_WL:
do i = 0, nct-1
SMW_CT(ct,i) = mw_sctran (SMW_MW(smw,i), system1,
"physical", axbits)
SMW_CTL(ct) = mw_sctran (SMW_MW(smw,0), "physical", system2,
axbits)
case SMW_PW, SMW_WP:
do i = 0, nct-1
SMW_CT(ct,i) = mw_sctran (SMW_MW(smw,i), system1,
system2, axbits)
}

return (ct)
end
#@(#) smwsmw.x 1.1
include "smw.h"


# SMW_SMW -- Set MCWS pointer

procedure smw_smw (smw, line, mw)

pointer smw #I SMW pointer
int line #I Physical line
pointer mw #I MWCS pointer

begin
if (SMW_NMW(smw) == 1)
SMW_MW(smw,0) = mw

else {
if (line < 1 || line > SMW_NSPEC(smw))
call error (1, "smw_smw: aperture not found")
SMW_MW(smw,(line-1)/SMW_NSPLIT) = mw
}
end
#@(#) smwswattrs.x 1.1
include
include "smw.h"


# SMW_SWATTRS -- Set spectrum attribute parameters
# This routine has the feature that if the coordinate system of a single
# spectrum in an EQUISPEC WCS is changed then the image WCS is changed
# to a MULTISPEC WCS.

procedure smw_swattrs (smw, index1, index2, ap, beam, dtype, w1, dw, nw, z,
aplow, aphigh, coeff)

pointer smw # SMW pointer
int index1 # Spectrum index
int index2 # Spectrum index
int ap # Aperture number
int beam # Beam number
int dtype # Dispersion type
double w1 # Starting coordinate
double dw # Coordinate interval
int nw # Number of valid pixels
double z # Redshift factor
real aplow[2], aphigh[2] # Aperture limits
char coeff[ARB] # Nonlinear coeff string

bool fp_equald()
int i, j, sz_val, strlen()
double a, b
pointer sp, str, val, mw
errchk smw_mw

define start_ 10

begin

call smark (sp)
call salloc (str, SZ_LINE, TY_CHAR)

start_
switch (SMW_FORMAT(smw)) {
case SMW_ND:
if (!IS_INDEFI(SMW_DTYPE(smw)) && (!fp_equald(w1,SMW_W1(smw)) ||
!fp_equald(dw,SMW_DW(smw)) || !fp_equald(z,SMW_Z(smw)))) {
call malloc (val, 15, TY_DOUBLE)
mw = SMW_MW(smw,0)
i = SMW_PDIM(smw)
j = SMW_PAXIS(smw,1)
call mw_gwtermd (mw, Memd[val], Memd[val+i], Memd[val+2*i], i)
Memd[val+j-1] = 1.
switch (dtype) {
case DCNO, DCLINEAR:
a = w1 / (1 + z)
b = dw / (1 + z)
case DCLOG:
a = log10 (w1 / (1 + z))
b = log10 ((w1 + (nw - 1) * dw) / w1) / (nw - 1)
case DCFUNC:
call error (1,
"Nonlinear functions not allowed for NSPEC format")
}
Memd[val+i+j-1] = a
Memd[val+2*i+(i+1)*(j-1)] = b
call mw_swtermd (mw, Memd[val], Memd[val+i], Memd[val+2*i], i)
call mfree (val, TY_DOUBLE)
}
SMW_DTYPE(smw) = dtype
SMW_NW(smw) = nw
SMW_W1(smw) = w1
SMW_DW(smw) = dw
SMW_Z(smw) = z

case SMW_ES:
# Check for any changes to the dispersion system.
if (dtype == DCFUNC) {
call smw_esms(smw)
goto start_
}
if (!IS_INDEFI(SMW_DTYPE(smw)) && (dtype != SMW_DTYPE(smw) ||
nw != SMW_NW(smw) || !fp_equald(w1,SMW_W1(smw)) ||
!fp_equald(dw,SMW_DW(smw)) || !fp_equald(z,SMW_Z(smw)))) {
if (SMW_NSPEC(smw) > 1 && index1 > 0) {
call smw_esms(smw)
goto start_
}
if (!fp_equald(w1,SMW_W1(smw)) || !fp_equald(dw,SMW_DW(smw)) ||
!fp_equald(z,SMW_Z(smw))) {
call malloc (val, 15, TY_DOUBLE)
mw = SMW_MW(smw,0)
i = SMW_PDIM(smw)
j = SMW_PAXIS(smw,1)
call mw_gwtermd (mw, Memd[val], Memd[val+i],
Memd[val+2*i], i)
Memd[val+j-1] = 1.
switch (dtype) {
case DCNO, DCLINEAR:
a = w1 / (1 + z)
b = dw / (1 + z)
case DCLOG:
a = log10 (w1 / (1 + z))
b = log10 ((w1 + (nw - 1) * dw) / w1) / (nw - 1)
}
Memd[val+i+j-1] = a
Memd[val+2*i+(i+1)*(j-1)] = b
call mw_swtermd (mw, Memd[val], Memd[val+i],
Memd[val+2*i], i)
call mfree (val, TY_DOUBLE)
}
}

SMW_DTYPE(smw) = dtype
SMW_NW(smw) = nw
SMW_W1(smw) = w1
SMW_DW(smw) = dw
SMW_Z(smw) = z

if (index1 > 0) {
Memi[SMW_APS(smw)+index1-1] = ap
Memi[SMW_BEAMS(smw)+index1-1] = beam
Memr[SMW_APLOW(smw)+2*index1-2] = aplow[1]
Memr[SMW_APHIGH(smw)+2*index1-2] = aphigh[1]
Memr[SMW_APLOW(smw)+2*index1-1] = aplow[2]
Memr[SMW_APHIGH(smw)+2*index1-1] = aphigh[2]
}

case SMW_MS:
# We can't use SPRINTF for the whole string because it can only
# handle a limited length and trucates long coefficient strings.
# Use STRCAT instead.

call smw_mw (smw, index1, index2, mw, i, j)
sz_val = strlen (coeff) + SZ_LINE
call salloc (val, sz_val, TY_CHAR)
call sprintf (Memc[str], SZ_LINE, "spec%d")
call pargi (i)
call sprintf (Memc[val], sz_val,
"%d %d %d %.14g %.14g %d %.14g %.2f %.2f")
call pargi (ap)
call pargi (beam)
call pargi (dtype)
if (dtype == DCLOG) {
call pargd (log10 (w1))
call pargd (log10 ((w1+(nw-1)*dw)/w1)/(nw-1))
} else {
call pargd (w1)
call pargd (dw)
}
call pargi (nw)
call pargd (z)
call pargr (aplow[1])
call pargr (aphigh[1])
if (coeff[1] != EOS) {
call strcat (" ", Memc[val], sz_val)
call strcat (coeff, Memc[val], sz_val)
}
call mw_swattrs (mw, 2, Memc[str], Memc[val])

if (SMW_APS(smw) != NULL)
Memi[SMW_APS(smw)+index1-1] = ap
}

call sfree (sp)
end
#@(#) t_specshift.x 1.1
include
include "smw.h"

# Function types.
define CHEBYSHEV 1 # CURFIT Chebyshev polynomial
define LEGENDRE 2 # CURFIT Legendre polynomial
define SPLINE3 3 # CURFIT cubic spline
define SPLINE1 4 # CURFIT linear spline
define PIXEL 5 # pixel coordinate array
define SAMPLE 6 # sampled coordinates


# T_SSHIFT -- Shift the spectral coordinates

procedure t_sshift ()

int list # Input list of spectra
double shift # Shift to apply
pointer aps # Aperture list
bool verbose # Verbose?

int ap, beam, dtype, nw
double w1, dw, z
real aplow[2], aphigh[2]
pointer sp, image, coeff, tmp, im, mw

bool clgetb()
double clgetd()
int imtopenp(), imtgetim()
pointer rng_open(), immap(), smw_openim()
errchk immap, smw_openim, smw_gwattrs, smw_swattrs, sshift

begin
call smark (sp)
call salloc (image, SZ_FNAME, TY_CHAR)
coeff = NULL

list = imtopenp ("spectra")
shift = clgetd ("shift")
call clgstr ("apertures", Memc[image], SZ_FNAME)
verbose = clgetb ("verbose")

iferr (aps = rng_open (Memc[image], INDEF, INDEF, INDEF))
call error (0, "Bad aperture list")

while (imtgetim (list, Memc[image], SZ_FNAME) != EOF) {
im = NULL
mw = NULL
iferr {
tmp = immap (Memc[image], READ_WRITE, 0); im = tmp
tmp = smw_openim (im); mw = tmp

switch (SMW_FORMAT(mw)) {
case SMW_ND:
call smw_gwattrs (mw, 1, 1, ap, beam, dtype,
w1, dw, nw, z, aplow, aphigh, coeff)
w1 = w1 + shift
call smw_swattrs (mw, 1, 1, ap, beam, dtype,
w1, dw, nw, z, aplow, aphigh, Memc[coeff])
if (verbose) {
call printf ("%s: shift = %g, %g --> %g\n")
call pargstr (Memc[image])
call pargd (shift)
call pargd (w1 - shift)
call pargd (w1)
}
case SMW_ES, SMW_MS:
call sshift (im, mw, Memc[image], aps, shift,
verbose)
}
} then
call erract (EA_WARN)

if (mw != NULL) {
call smw_saveim (mw, im)
call smw_close (mw)
}
if (im != NULL)
call imunmap (im)
}

call rng_close (aps)
call imtclose (list)
call mfree (coeff, TY_CHAR)
call sfree (sp)
end


# SSHIFT -- Shift coordinate zero point of selected aperture in
# MULTISPEC and EQUISPEC images.

procedure sshift (im, mw, image, aps, shift, verbose)

pointer im # IMIO pointer
pointer mw # MWCS pointer
char image[ARB] # Image name
pointer aps # Aperture range list
double shift # Shift to add
bool verbose # Verbose?

int i, ap, beam, dtype, nw, naps
double w1, dw, z
real aplow[2], aphigh[2]
pointer coeff, coeffs
bool rng_elementi()
errchk sshift1

begin
coeff = NULL
coeffs = NULL

# Go through each spectrum and change the selected apertures.
naps = 0
do i = 1, SMW_NSPEC(mw) {
# Get aperture info
iferr (call smw_gwattrs (mw, i, 1, ap, beam, dtype, w1, dw, nw, z,
aplow, aphigh, coeff))
break

# Check if aperture is to be changed
if (!rng_elementi (aps, ap))
next

# Apply shift
w1 = w1 + shift
if (dtype == 2)
call sshift1 (shift, coeff)

call smw_swattrs (mw, i, 1, ap, beam, dtype, w1, dw, nw, z,
aplow, aphigh, Memc[coeff])

# Make record
if (verbose) {
if (naps == 1) {
call printf ("%s: shift = %g\n")
call pargstr (image)
call pargd (shift)
}
call printf (" Aperture %d: %g --> %g\n")
call pargi (ap)
call pargd (w1 - shift)
call pargd (w1)
}
}

call mfree (coeff, TY_CHAR)
call mfree (coeffs, TY_DOUBLE)
end


# SSHIFT1 -- Shift coordinate zero point of nonlinear functions.

procedure sshift1 (shift, coeff)

double shift # Shift to add
pointer coeff # Attribute function coefficients

int i, j, ip, nalloc, ncoeff, type, order, fd
double dval
pointer coeffs
int ctod(), stropen()
errchk stropen

begin
if (coeff == NULL)
return
if (Memc[coeff] == EOS)
return

coeffs = NULL
ncoeff = 0
ip = 1
while (ctod (Memc[coeff], ip, dval) > 0) {
if (coeffs == NULL) {
nalloc = 10
call malloc (coeffs, nalloc, TY_DOUBLE)
} else if (ncoeff == nalloc) {
nalloc = nalloc + 10
call realloc (coeffs, nalloc, TY_DOUBLE)
}
Memd[coeffs+ncoeff] = dval
ncoeff = ncoeff + 1
}
ip = ip + SZ_LINE
call realloc (coeff, ip, TY_CHAR)
call aclrc (Memc[coeff], ip)
fd = stropen (Memc[coeff], ip, NEW_FILE)

ip = 0
while (ip < ncoeff) {
if (ip > 0)
call fprintf (fd, " ")
Memd[coeffs+ip+1] = Memd[coeffs+ip+1] + shift
type = nint (Memd[coeffs+ip+2])
order = nint (Memd[coeffs+ip+3])
call fprintf (fd, "%.3g %g %d %d")
call pargd (Memd[coeffs+ip])
call pargd (Memd[coeffs+ip+1])
call pargi (type)
call pargi (order)
switch (type) {
case CHEBYSHEV, LEGENDRE:
j = 6 + order
case SPLINE3:
j = 9 + order
case SPLINE1:
j = 7 + order
case PIXEL:
j = 4 + order
case SAMPLE:
j = 5 + order
}
do i = 4, j-1 {
call fprintf (fd, " %g")
call pargd (Memd[coeffs+ip+i])
}
ip = ip + j
}
call strclose (fd)

call mfree (coeffs, TY_DOUBLE)
end
#@(#) units.x 1.1
include
include "units.h"


# UN_OPEN -- Open units package
# It is allowed to open an unknown unit type

pointer procedure un_open (units)

char units[ARB] # Units string
pointer un # Units pointer returned

begin
call calloc (un, UN_LEN, TY_STRUCT)
iferr (call un_decode (un, units)) {
call un_close (un)
call erract (EA_ERROR)
}
return (un)
end


# UN_CLOSE -- Close units package

procedure un_close (un)

pointer un # Units pointer

begin
call mfree (un, TY_STRUCT)
end


# UN_COPY -- Copy units pointer

procedure un_copy (un1, un2)

pointer un1, un2 # Units pointers

begin
if (un2 == NULL)
call malloc (un2, UN_LEN, TY_STRUCT)
call amovi (Memi[un1], Memi[un2], UN_LEN)
end


# UN_DECODE -- Decode units string and set up units structure.
# The main work is done in UN_DECODE1 so that the units string may
# be recursive; i.e. the units string may contain other units strings.
# In particular, this is required for the velocity units to specify
# a reference wavelength.

procedure un_decode (un, units)

pointer un # Units pointer
char units[ARB] # Units string

bool streq()
pointer sp, units1, temp, un1, un2
errchk un_decode1, un_ctranr

begin
if (streq (units, UN_USER(un)))
return

call smark (sp)
call salloc (units1, SZ_LINE, TY_CHAR)
call salloc (temp, UN_LEN, TY_STRUCT)

# Save a copy to restore in case of an error.
call un_copy (un, temp)

iferr {
# Decode the primary units
call un_decode1 (un, units, Memc[units1], SZ_LINE)

# Decode velocity reference wavelength if necessary.
if (UN_CLASS(un) == UN_VEL || UN_CLASS(un) == UN_DOP) {
call salloc (un1, UN_LEN, TY_STRUCT)
call un_decode1 (un1, Memc[units1], Memc[units1], SZ_LINE)
if (UN_CLASS(un1) == UN_VEL || UN_CLASS(un1) == UN_DOP)
call error (1,
"Velocity reference units may not be velocity")
call salloc (un2, UN_LEN, TY_STRUCT)
call un_decode1 (un2, "angstroms", Memc[units1], SZ_LINE)
call un_ctranr (un1, un2, UN_VREF(un), UN_VREF(un), 1)
}
} then {
call un_copy (temp, un)
call sfree (sp)
call erract (EA_ERROR)
}

call sfree (sp)
end


# UN_DECODE1 -- Decode units string and set up units structure.
# Return any secondary units string. Unknown unit strings are allowed.

procedure un_decode1 (un, units, units1, sz_units1)

pointer un # Units pointer
char units[ARB] # Units string
char units1[sz_units1] # Secondary units string to return
int sz_units1 # Size of secondary units string

int unlog, uninv, untype
int i, j, nscan(), strdic()
pointer sp, str

int class[UN_NUNITS]
real scale[UN_NUNITS]
data class /UN_WAVE,UN_WAVE,UN_WAVE,UN_WAVE,UN_WAVE,UN_WAVE,UN_WAVE,
UN_FREQ,UN_FREQ,UN_FREQ,UN_FREQ,UN_VEL,UN_VEL,
UN_ENERGY,UN_ENERGY,UN_ENERGY,UN_DOP/
data scale /UN_ANG,UN_NM,UN_MMIC,UN_MIC,UN_MM,UN_CM,UN_M,UN_HZ,UN_KHZ,
UN_MHZ,UN_GHZ,UN_MPS,UN_KPS,UN_EV,UN_KEV,UN_MEV,UN_Z/

begin
call smark (sp)
call salloc (str, SZ_FNAME, TY_CHAR)

call strcpy (units, Memc[str], SZ_FNAME)
call strlwr (Memc[str])
call sscan (Memc[str])
untype = 0
unlog = NO
uninv = NO
do i = 1, 3 {
call gargwrd (Memc[str], SZ_FNAME)
if (nscan() != i)
break

j = strdic (Memc[str], Memc[str], SZ_FNAME, UN_DIC)
if (j > UN_NUNITS) {
j = j - UN_NUNITS
if (j == 1) {
if (unlog == YES)
break
unlog = YES
} else if (j == 2) {
if (uninv == YES)
break
uninv = YES
}
} else {
if (class[j] == UN_VEL || class[j] == UN_DOP) {
call gargr (UN_VREF(un))
call gargstr (units1, sz_units1)
if (nscan() != i+2)
call error (1, "Error in velocity reference wavelength")
} else
UN_VREF(un) = 0.
untype = j
break
}
}

if (untype == 0) {
UN_TYPE(un) = 0
UN_CLASS(un) = UN_UNKNOWN
UN_LABEL(un) = EOS
call strcpy (units, UN_UNITS(un), SZ_UNITS)
} else {
UN_TYPE(un) = untype
UN_CLASS(un) = class[untype]
UN_LOG(un) = unlog
UN_INV(un) = uninv
UN_SCALE(un) = scale[untype]
UN_LABEL(un) = EOS
UN_UNITS(un) = EOS
call strcpy (units, UN_USER(un), SZ_UNITS)

if (unlog == YES)
call strcat ("Log ", UN_LABEL(un), SZ_UNITS)
if (uninv == YES)
call strcat ("inverse ", UN_UNITS(un), SZ_UNITS)
call strcat (Memc[str], UN_UNITS(un), SZ_UNITS)
switch (class[j]) {
case UN_WAVE:
if (uninv == NO)
call strcat ("Wavelength", UN_LABEL(un), SZ_UNITS)
else
call strcat ("Wavenumber", UN_LABEL(un), SZ_UNITS)
case UN_FREQ:
call strcat ("Frequency", UN_LABEL(un), SZ_UNITS)
case UN_VEL:
call strcat ("Velocity", UN_LABEL(un), SZ_UNITS)
case UN_ENERGY:
call strcat ("Energy", UN_LABEL(un), SZ_UNITS)
case UN_DOP:
call strcat ("Redshift", UN_LABEL(un), SZ_UNITS)
}
}

call sfree (sp)
end


# UN_COMPARE -- Compare two units

bool procedure un_compare (un1, un2)

pointer un1, un2 # Units pointers to compare
bool strne()

begin
if (strne (UN_UNITS(un1), UN_UNITS(un2)))
return (false)
if (strne (UN_LABEL(un1), UN_LABEL(un2)))
return (false)
if (UN_VREF(un1) != UN_VREF(un2))
return (false)
return (true)
end


# UN_CTRANR -- Transform units
# Error is returned if the transform cannot be made

procedure un_ctranr (un1, un2, val1, val2, nvals)

pointer un1 # Input units pointer
pointer un2 # Output units pointer
real val1[nvals] # Input values
real val2[nvals] # Output values
int nvals # Number of values

int i
real s, v, z
bool un_compare()

begin
if (un_compare (un1, un2)) {
call amovr (val1, val2, nvals)
return
}

if (UN_CLASS(un1) == UN_UNKNOWN || UN_CLASS(un2) == UN_UNKNOWN)
call error (1, "Cannot convert between selected units")

call amovr (val1, val2, nvals)

s = UN_SCALE(un1)
if (UN_LOG(un1) == YES)
do i = 1, nvals
val2[i] = 10. ** val2[i]
if (UN_INV(un1) == YES)
do i = 1, nvals
val2[i] = 1. / val2[i]
switch (UN_CLASS(un1)) {
case UN_WAVE:
do i = 1, nvals
val2[i] = val2[i] / s
case UN_FREQ:
do i = 1, nvals
val2[i] = s / val2[i]
case UN_VEL:
v = UN_VREF(un1)
do i = 1, nvals {
z = val2[i] / s
val2[i] = sqrt ((1 + z) / (1 - z)) * v
}
case UN_ENERGY:
do i = 1, nvals
val2[i] = s / val2[i]
case UN_DOP:
v = UN_VREF(un1)
do i = 1, nvals
val2[i] = (val2[i] / s + 1) * v
}

s = UN_SCALE(un2)
switch (UN_CLASS(un2)) {
case UN_WAVE:
do i = 1, nvals
val2[i] = val2[i] * s
case UN_FREQ:
do i = 1, nvals
val2[i] = s / val2[i]
case UN_VEL:
v = UN_VREF(un2)
do i = 1, nvals {
z = (val2[i] / v) ** 2
val2[i] = (z - 1) / (z + 1) * s
}
case UN_ENERGY:
do i = 1, nvals
val2[i] = s / val2[i]
case UN_DOP:
v = UN_VREF(un2)
do i = 1, nvals
val2[i] = (val2[i] / v - 1) * s
}
if (UN_INV(un2) == YES)
do i = 1, nvals
val2[i] = 1. / val2[i]
if (UN_LOG(un2) == YES)
do i = 1, nvals
val2[i] = log10 (val2[i])
end


# UN_CHANGER -- Change units
# Error is returned if the conversion cannot be made

procedure un_changer (un, units, vals, nvals, update)

pointer un # Units pointer (may be changed)
char units[ARB] # Desired units
real vals[nvals] # Values
int nvals # Number of values
int update # Update units pointer?

bool streq(), un_compare()
pointer un1, un_open()
errchk un_open, un_ctranr

begin

# Check for same unit string
if (streq (units, UN_USER(un)))
return

# Check for error in units string, or the same units.
un1 = un_open (units)
if (un_compare (un1, un)) {
call strcpy (units, UN_USER(un), SZ_UNITS)
call un_close (un1)
return
}

iferr {
call un_ctranr (un, un1, vals, vals, nvals)
if (update == YES)
call un_copy (un1, un)
call un_close(un1)
} then {
call un_close(un1)
call erract (EA_ERROR)
}
end


# UN_CTRAND -- Transform units
# Error is returned if the transform cannot be made

procedure un_ctrand (un1, un2, val1, val2, nvals)

pointer un1 # Input units pointer
pointer un2 # Output units pointer
double val1[nvals] # Input values
double val2[nvals] # Output values
int nvals # Number of values

int i
double s, v, z
bool un_compare()

begin
if (un_compare (un1, un2)) {
call amovd (val1, val2, nvals)
return
}

if (UN_CLASS(un1) == UN_UNKNOWN || UN_CLASS(un2) == UN_UNKNOWN)
call error (1, "Cannot convert between selected units")

call amovd (val1, val2, nvals)

s = UN_SCALE(un1)
if (UN_LOG(un1) == YES)
do i = 1, nvals
val2[i] = 10. ** val2[i]
if (UN_INV(un1) == YES)
do i = 1, nvals
val2[i] = 1. / val2[i]
switch (UN_CLASS(un1)) {
case UN_WAVE:
do i = 1, nvals
val2[i] = val2[i] / s
case UN_FREQ:
do i = 1, nvals
val2[i] = s / val2[i]
case UN_VEL:
v = UN_VREF(un1)
do i = 1, nvals {
z = val2[i] / s
val2[i] = sqrt ((1 + z) / (1 - z)) * v
}
case UN_ENERGY:
do i = 1, nvals
val2[i] = s / val2[i]
case UN_DOP:
v = UN_VREF(un1)
do i = 1, nvals
val2[i] = (val2[i] / s + 1) * v
}

s = UN_SCALE(un2)
switch (UN_CLASS(un2)) {
case UN_WAVE:
do i = 1, nvals
val2[i] = val2[i] * s
case UN_FREQ:
do i = 1, nvals
val2[i] = s / val2[i]
case UN_VEL:
v = UN_VREF(un2)
do i = 1, nvals {
z = (val2[i] / v) ** 2
val2[i] = (z - 1) / (z + 1) * s
}
case UN_ENERGY:
do i = 1, nvals
val2[i] = s / val2[i]
case UN_DOP:
v = UN_VREF(un2)
do i = 1, nvals
val2[i] = (val2[i] / v - 1) * s
}
if (UN_INV(un2) == YES)
do i = 1, nvals
val2[i] = 1. / val2[i]
if (UN_LOG(un2) == YES)
do i = 1, nvals
val2[i] = log10 (val2[i])
end


# UN_CHANGED -- Change units
# Error is returned if the conversion cannot be made

procedure un_changed (un, units, vals, nvals, update)

pointer un # Units pointer (may be changed)
char units[ARB] # Desired units
double vals[nvals] # Values
int nvals # Number of values
int update # Update units pointer?

bool streq(), un_compare()
pointer un1, un_open()
errchk un_open, un_ctrand

begin

# Check for same unit string
if (streq (units, UN_USER(un)))
return

# Check for error in units string, or the same units.
un1 = un_open (units)
if (un_compare (un1, un)) {
call strcpy (units, UN_USER(un), SZ_UNITS)
call un_close (un1)
return
}

iferr {
call un_ctrand (un, un1, vals, vals, nvals)
if (update == YES)
call un_copy (un1, un)
call un_close(un1)
} then {
call un_close(un1)
call erract (EA_ERROR)
}
end
#@(#) usercoord.x 1.1
include
include "smw.h"
include "units.h"

# USERCOORD -- Set user coordinates

procedure usercoord (sh, key, w1, u1, w2, u2)

pointer sh
int key
double w1, u1, w2, u2

int i, format, ap, beam, dtype, nw
double shift, wa, wb, ua, ub, w0, dw, z, smw_c1trand()
real aplow[2], aphigh[2]
pointer coeff, smw, mw, ct, smw_sctran()
errchk smw_sctran

begin
coeff = NULL
smw = MW(sh)
mw = SMW_MW(smw,0)
format = SMW_FORMAT(smw)

iferr {
call un_ctrand (UN(sh), MWUN(sh), w1, wa, 1)
call un_ctrand (UN(sh), MWUN(sh), u1, ua, 1)

call smw_gwattrs (MW(sh), APINDEX(sh), LINDEX(sh,2),
ap, beam, dtype, w0, dw, nw, z, aplow, aphigh, coeff)

switch (key) {
case 'd':
wa = wa * (1 + z)
switch (UN_CLASS(MWUN(sh))) {
case UN_WAVE:
z = (wa - ua) / ua
case UN_FREQ, UN_ENERGY:
z = (ua - wa) / wa
default:
call error (1, "Inappropriate coordinate units")
}
case 'z':
shift = ua - wa
w0 = w0 + shift
if (dtype == 2)
call sshift1 (shift, coeff)
case 'l':
call un_ctrand (UN(sh), MWUN(sh), w2, wb, 1)
call un_ctrand (UN(sh), MWUN(sh), u2, ub, 1)

switch (format) {
case SMW_ND:
i = 2 ** (SMW_PAXIS(smw,1) - 1)
ct = smw_sctran (smw, "world", "physical", i)
if (dtype == DCLOG) {
wa = log10 (wa)
wb = log10 (wb)
}
wa = smw_c1trand (ct, wa)
wb = smw_c1trand (ct, wb)
case SMW_ES, SMW_MS:
ct = smw_sctran (smw, "world", "physical", 3)
if (dtype == DCLOG) {
wa = log10 (wa)
wb = log10 (wb)
}
call smw_c2trand (ct, wa, double (ap), wa, shift)
call smw_c2trand (ct, wb, double (ap), wb, shift)
}
call smw_ctfree (ct)

dw = (ub - ua) / (wb - wa)
w0 = ua - (wa - 1) * dw
dtype = 0
if (UNITS(sh) == EOS) {
call mw_swattrs (mw, SMW_PAXIS(smw,1),
"label", "Wavelength")
call mw_swattrs (mw, SMW_PAXIS(smw,1),
"units", "angstroms")
}
default:
call error (1, "Unknown correction")
}

call smw_swattrs (smw, LINDEX(sh,1), 1, ap, beam, dtype, w0,
dw, nw, z, aplow, aphigh, Memc[coeff])
if (smw != MW(sh)) {
CTLW1(sh) = NULL
CTWL1(sh) = NULL
MW(sh) = smw
}

DC(sh) = dtype
call shdr_system (sh, "world")
if (UN_CLASS(UN(sh)) == UN_UNKNOWN)
call un_copy (MWUN(sh), UN(sh))
} then
call erract (EA_WARN)

call mfree (coeff, TY_CHAR)
end



CMP_DELETE (24Jan96) source CMP_DELETE (24Jan96)



NAME
cmp_delete -- Delete a component


DESCRIPTION



ROUTINE DETAILS

procedure cmp_delete (cmp_list, idstr)


ARGUMENTS
cmp_list (input: pointer)
Component list

idstr[ARB] (input: char)
Id of component to delete



CMP_GET (24Jan96) source CMP_GET (24Jan96)



NAME
cmp_get -- Get a component based on an id string


DESCRIPTION



ROUTINE DETAILS

pointer procedure cmp_get (cmp_list, idstr, type)


ARGUMENTS
cmp_list (input: pointer)
Component list

idstr[ARB] (input: char)
id

type[ARB] (input: char)
type of component

RETURNS



CMP_INQUIRE (24Jan96) source CMP_INQUIRE (24Jan96)



NAME
cmp_inquire -- Get information about a component


DESCRIPTION



ROUTINE DETAILS

procedure cmp_inquire (cmp_list, idstr)


ARGUMENTS
cmp_list (input: pointer)
component list

idstr[ARB] (input: char)
Id of component to find information



CMP_MARK (24Jan96) source CMP_MARK (24Jan96)



NAME
cmp_mark -- Send component values and create component marker


DESCRIPTION
This will need to change when component drawing is finalized. The
important point to note here is that sending the parameters of a
component should always be accompanied by a marking of the
component. The two should be considered an "atomic" operation.


ROUTINE DETAILS

procedure cmp_mark (as, id)


ARGUMENTS
as (input: pointer)
Aspect object

id (input: int)
Id of component to mark



CMP_SET (24Jan96) source CMP_SET (24Jan96)



NAME
cmp_set -- Set parameters of a component


DESCRIPTION



ROUTINE DETAILS

procedure cmp_set_from_params (cmp_list, idstr, type, params)


ARGUMENTS
cmp_list (input: pointer)
Component list

idstr[ARB] (input: char)
Id of component

type[ARB] (input: char)
Type of component

params[ARB] (input: char)
String containing parameter values



CMP_SET_FROM_BOX (24Jan96) source CMP_SET_FROM_BOX (24Jan96)



NAME
cmp_set_from_box -- Set component parameters from a box


DESCRIPTION



ROUTINE DETAILS

procedure cmp_set_from_box (as, idstr, type, left, right, bottom, top,
source)

ARGUMENTS
as (input: pointer)
Aspect Object

idstr[ARB] (input: char)
Id of component to set

right (input: double)
orizontal range of box in NDC

left (input: double)
orizontal range of box in NDC

bottom (input: double)
Vertical range of box in NDC

top (input: double)
Vertical range of box in NDC

source[ARB] (input: char0
Source of NDC coordinates



COM_CLIENT (24Jan96) source COM_CLIENT (24Jan96)



NAME
com_client -- Execute client-related activities


DESCRIPTION



ROUTINE DETAILS

procedure com_client (as, command)


ARGUMENTS
as (input: pointer)
the ASpec object

command[ARB] (input: char)
The command to execute



COM_CMP (24Jan96) source COM_CMP (24Jan96)



NAME
com_cmp -- Execute client-related activities


DESCRIPTION



ROUTINE DETAILS

procedure com_cmp (as, command)


ARGUMENTS
as (input: pointer)
The ASpect object

command[ARB] (input: int)
Command to execute



COM_DB (24Jan96) source COM_DB (24Jan96)



NAME
com_db -- Execute component database-related activities


DESCRIPTION



ROUTINE DETAILS

procedure com_db (as, command)


ARGUMENTS
as (input: pointer)
the ASpec object

command[ARB] (input: char)
command to execute



COM_FILE (24Jan96) source COM_FILE (24Jan96)



NAME
com_file -- Execute file-related activities


DESCRIPTION



ROUTINE DETAILS

procedure com_file (as, command)


ARGUMENTS
as (input: pointer)
The ASpect object

command[ARB] (input: char)
command to execute



COM_FIT (24Jan96) source COM_FIT (24Jan96)



NAME
com_fit -- Execute fit-related activities


DESCRIPTION



ROUTINE DETAILS

procedure com_fit (as, command)


ARGUMENTS
as (input: pointer)
The ASpect object

command[ARB] (input: char)
Command to execute



COM_GTERM (24Jan96) source COM_GTERM (24Jan96)



NAME
com_gterm -- Execute gterm-related activities


DESCRIPTION



ROUTINE DETAILS

procedure com_gterm (as, command)


ARGUMENTS
as (input: pointer)
the ASpec object

command[ARB] (input: char)
command to execute



COM_RANGE (24Jan96) source COM_RANGE (24Jan96)



NAME
com_range -- Execute fit range-related activities


DESCRIPTION



ROUTINE DETAILS

procedure com_range (as, command)


ARGUMENTS
as (input: pointer)
the ASpec object

command[ARB] (input: char)
Command string



COM_SP (224Jan96) source COM_SP (224Jan96)



NAME
com_sp -- Execute spectrum file-related activities


DESCRIPTION



ROUTINE DETAILS

procedure com_sp (as, command)


ARGUMENTS
as (input: pointer)
the ASpec object

command[ARB] (input: char)
command to execute



FIT_EXECUTE (24Jan96) source/tasks/vuespec FIT_EXECUTE (24Jan96)



NAME
fit_execute -- Perform the fit operation


DESCRIPTION


ROUTINE DETAILS

procedure fit_execute (as)
ARGUMENTS
as (input: pointer)
the ASpec object



FIT_SET (6Feb96) source/tasks/vuespec FIT_SET (6Feb96)



NAME
fit_set -- Set the fitting parameters


DESCRIPTION



ROUTINE DETAILS

procedure fit_set (ft, p)
ARGUMENTS
ft (input: pointer)
fit object

p (input: pointer)
params object



PARSE (24Jan96) source PARSE (24Jan96)



NAME
parse -- Parse a colon command


DESCRIPTION



ROUTINE DETAILS

procedure parse (as, command)


ARGUMENTS
as (input: pointer)
the ASpec object

command[ARB] (input: char)
command line



PLOT_SPEC (24Jan96) source PLOT_SPEC (24Jan96)



NAME
plot_spec -- Plot the spectra


DESCRIPTION



ROUTINE DETAILS

procedure plot_spec (as)


ARGUMENTS
as (input: pointer)
ASpec object



RANGE_ADD (24Jan96) source RANGE_ADD (24Jan96)



NAME
range_add -- Add a domain range to the fit parameters


DESCRIPTION



ROUTINE DETAILS

procedure range_add (as, source, x1, x2)


ARGUMENTS
as (input: pointer)
ASpec object

source[ARB] (input: char)
Name of the source plot for coords

x1 (input: double)
Domain range in NDC (not ordered)

x2 (input: double)
Domain range in NDC (not ordered)



RANGE_DEL (24Jan96) source RANGE_DEL (24Jan96)



NAME
range_del -- Delete a domain range


DESCRIPTION



ROUTINE DETAILS

procedure range_del (as, plot_name, x)


ARGUMENTS
as (input: pointer)
ASpec object

plot_name[ARB] (input: char)
Plot source for the coordinat

x (input: double)
NDC Coordinate within range to delete



T_ASPECT (24Jan96) source T_ASPECT (24Jan96)



NAME
aspect -- Spectrum fitting


DESCRIPTION
interactive spectrum modelling task. Initializes "aspect" object &
calls event loop.


ROUTINE DETAILS

procedure t_aspect()


ARGUMENTS



UNDO_SETS (24Jan96) source UNDO_SETS (24Jan96)



NAME
un_cmp_list -- Setup to restore the component list


DESCRIPTION



ROUTINE DETAILS

procedure un_cmp_list (un, c)


ARGUMENTS
un (input: pointer)
UNDO ring

c (input: pointer)
Component list



UNDO (24Jan96) source UNDO (24Jan96)



NAME
undo -- Undo an operation


DESCRIPTION



ROUTINE DETAILS

procedure undo (as)


ARGUMENTS
as (input: pointer)
the ASpec object



UNDO_STRUCT (24Jan96) source UNDO_STRUCT (24Jan96)



NAME
un_alloc -- Allocate the UNDO ring
un_free -- Free the ring
unop_alloc -- Allocate the UNDO operation.
unop_free -- Free the operation
un_push -- Push an operation onto the ring.
un_pop -- Pop an operation from the ring.
un_get -- Get operation on ring without removing it


DESCRIPTION



ROUTINE DETAILS

pointer procedure un_alloc (size)


ARGUMENTS
size (input: int)
size of undo ring


procedure un_free (r)


ARGUMENTS
r (input/output: pointer)
Undo ring, NULL on exit


pointer procedure unop_alloc ()


ARGUMENTS

RETURNS
UNDO operation object


pointer procedure un_pop (r)


ARGUMENTS
r (input: pointer)
UNDO ring


pointer procedure un_get (r)


ARGUMENTS
r (input: pointer)
UNDO ring

RETURNS
pointer




GIO_WCS (25Jan96) source GIO_WCS (25Jan96)



NAME
gio_wcs -- Manage multiple GIO WCS contexts
ui_alloc_wcs -- Create a WCS object
ui_get_wcs -- Return the current WCS object
ui_set_wcs -- Set the current WCS using a specified WCS object
ui_switch_wcs -- Switch between WCS objects
ui_free_wcs -- Destroy the specified WCS object


DESCRIPTION
====================WARNING==================== THE BELOW
REPRESENTS AN INTERFACE VIOLATION WITH RELATION TO THE IRAF GIO
INTERFACE. TAKE CAREFUL NOTE OF IF AND WHEN GIO CHANGES AFFECT HOW
IT STORES WCS INFORMATION. WHEN IT CHANGES, THE BELOW ROUTINES
WILL PROBABLY BREAK ===============================================

In order to handle multiple graphics windows available to GIO, when
GIO itself can only handle one at a time, the below routines save
and restore GIO world coordinate (WCS) information, so the
application can retain this knowledge while switching between
graphic windows.


USAGE
In general, first create the needed WCS objects:
wcs = ui_alloc_wcs (NULL)

Then, each time a switch is made between graphic windows,
retrieve the current wcs information and set the one desired:

call ui_get_wcs (gp, wcs_old)
call ui_set_wcs (gp, wcs_new)

or

call ui_switch_wcs (gp, wcs_new, wcs_old)

Finally, when the application is finished with the objects,
remove them.

call ui_free_wcs (wcs)

UI_ALLOC_WCS COMMENTS
It should be noted that ui_alloc_wcs can optionally fill the WCS
object with the current GIO WCS information if passed the
graphics descriptor. If the argument is NULL, an empty WCS
object is returned.


ROUTINE DETAILS

pointer procedure ui_alloc_wcs


ARGUMENTS
gp (input: pointer)
graphics descriptor

RETURNS
WCS object



procedure ui_get_wcs (gp, wcs)


ARGUMENTS
gp (input: pointer)
graphics descriptor

wcs (input: pointer)
WCS Object


procedure ui_set_wcs (gp, wcs)


ARGUMENTS
gp (input: pointer)
graphics descriptor

wcs (input: pointer)
WCS WCS state to return to


procedure ui_free_wcs (wcs)


ARGUMENTS
wcs (input/output: pointer)
WCS state, NULL on exit


procedure ui_switch_wcs (gp, wcs_new, wcs_old)


ARGUMENTS
gp (input: pointer)
graphics descriptor

wcs_new (input: pointer)
WCS information to set GIO

wcs_old (input: pointer)
Return the current WCS information


procedure ui_copy_wcs (input, output)


ARGUMENTS
input (input: pointer)
WCS to copy

output (input: pointer)
WCS to copy to



UI_CONVERT (25Jan96) source UI_CONVERT (25Jan96)



NAME
ui_toint -- Convert to internal units
ui_fromint -- Convert from internal units


DESCRIPTION



ROUTINE DETAILS

procedure ui_fromint (funit, dunit, fi, di, f, d)


ARGUMENTS
funit (input: int)
Flux units to convert to

dunit (input: int)
Dispersion units to convert to

fi (input: real)
Flux to convert

di (input: double)
Dispersion to convert

f (output: real)
Converted flux

d (output: double)
Converted dispersion


procedure ui_toint (funit, dunit, f, d, fi, di)


ARGUMENTS
funit (input: int)
Flux units to convert to

dunit (input: int)
Dispersion units to convert to

f (output: real)
Converted flux

d (output: double)
Converted dispersion

fi (input: real)
Internal Flux

di (input: double)
Internal Dispersion



UI_ALLOC (25Jan96) source UI_ALLOC (25Jan96)



NAME
ui_alloc -- Setup the UI parameters


DESCRIPTION



ROUTINE DETAILS

procedure ui_alloc (device, gui_file)


ARGUMENTS
device[ARB] (input: char)
Graphics device name

gui_file[ARB] (input: char)
GUI initialization file



UI_ALLOC_EVENT (25Jan96) source UI_ALLOC_EVENT (25Jan96)



NAME
ui_alloc_event -- Allocate a UI object


DESCRIPTION



ROUTINE DETAILS

pointer procedure ui_alloc_event ()


RETURNS
the object pointer



UI_BUSY (25Jan96) source UI_BUSY (25Jan96)



NAME
ui_busy -- Send busy status


DESCRIPTION



ROUTINE DETAILS

procedure ui_busy (state)


ARGUMENTS
sate (input: bool)
Busy?



UI_CHDIR (25Jan96) source UI_CHDIR (25Jan96)



NAME
ui_chdir -- Change the directory and report on subdirectories


DESCRIPTION



ROUTINE DETAILS

procedure ui_chdir (dir)


ARGUMENTS
dir[ARB] (input: char)
The directory to change to



UI_CLIP_VP (25Jan96) source UI_CLIP_VP (25Jan96)



NAME
ui_clip_vp -- Clip coordinates to viewport


DESCRIPTION



ROUTINE DETAILS

procedure ui_clip_vp (as, x, y, s, sys, cx, cy)


ARGUMENTS
as (input: pointer)
ASpec object

x (input: double)
Coordinates to clip

y (input: double)
Coordinates to clip

s[ARB] (input: char)
Source of coordinates

sys (input: int)
Clip in NDC or WORLD coordinates

cx (output: double)
Clipped coordinates

cy (output: double)
Clipped coordinates



UI_ERRMSG (25Jan96) source UI_ERRMSG (25Jan96)



NAME
ui_errmsg -- Send an error message to the gui


DESCRIPTION



ROUTINE DETAILS

procedure ui_errmsg (msg)


ARGUMENTS
msg[ARB] (input: char)
error message



UI_EVENT_LOOP (25Jan96) source UI_EVENT_LOOP (25Jan96)



NAME
ui_event_loop -- Handle events from the GUI


DESCRIPTION
This routine handles all events generated by the GUI.


ROUTINE DETAILS

procedure ui_event_loop (as)


ARGUMENTS
as (input: pointer)
ASpec object



UI_FATAL_ERROR (25Jan96) source UI_FATAL_ERROR (25Jan96)



NAME
ui_fatal_error -- Handle fatal errors


DESCRIPTION



ROUTINE DETAILS

procedure ui_fatal_error (err, next_handle)


ARGUMENTS
err (input: int)
error code

next_handle (output: int)
next handler to call



UI_FIT_MARK (25Jan96) source UI_FIT_MARK (25Jan96)



NAME
ui_fit_mark - Plot the current fit


DESCRIPTION



ROUTINE DETAILS

procedure ui_fit_mark (as)


ARGUMENTS
as (input: pointer)
ASpec object



UI_FREE_EVENT (25Jan96) source UI_FREE_EVENT (25Jan96)



NAME
ui_free_event -- Free a UI object


DESCRIPTION



ROUTINE DETAILS

procedure ui_free_event (o)


ARGUMENTS
o (input/output: pointer)
The object. NULL on return



UI_FROM_NDC (25Jan96) source UI_FROM_NDC (25Jan96)



NAME
ui_from_ndc -- Convert from NDC to User coordinates


DESCRIPTION



ROUTINE DETAILS

procedure ui_from_ndc (as, x, y, source, ux, uy)


ARGUMENTS
as (input: pointer)
ASpec object

x (input: double)
NDC coordinates

y (input: double)
NDC coordinates

source[ARB] (input: char)
Source of the NDC coordinates

ux (output: double)
User coordinates

uy (output: double)
User coordinates



UI_HIGHLIGH_CMP (25Jan96) source UI_HIGHLIGH_CMP (25Jan96)



NAME
ui_highlight_cmp -- Highlight a single component


DESCRIPTION
Stopgap routine until components are markers.


ROUTINE DETAILS

procedure ui_highlight_cmp (as)


ARGUMENTS
as (input: pointer)
ASpec object



UI_INIT (25Jan96) source UI_INIT (25Jan96)



NAME
ui_init -- Do some basic GUI initialization


DESCRIPTION
This routine initializes the UI common block, reads in relevant CL
parameters, and informs the gui of the current state of the task.

If the task has been re-run and the user wishes to preserve the
state from a previous run, simply not calling this routine will
accomplish this.


ROUTINE DETAILS

procedure ui_init (as)


ARGUMENTS
as (input: pointer)
ASpec onbject



UI_LS (25Jan96) source UI_LS (25Jan96)



NAME
ui_ls -- Get a directory listing


DESCRIPTION



ROUTINE DETAILS

procedure ui_ls (filter)


ARGUMENTS
filter[ARB] (input: char)
Filter for the directory listing



UI_MARK_RANGE (25Jan96) source UI_MARK_RANGE (25Jan96)



NAME
ui_mark_range -- Mark the specified range


DESCRIPTION



ROUTINE DETAILS

procedure ui_mark_range (di, df)


ARGUMENTS
di (input: double)
Dispersion (angstroms) range to mark

df (input: double)
Dispersion (angstroms) range to mark



UI_MATH_ERROR (25Jan96) source UI_MATH_ERROR (25Jan96)



NAME
ui_math_error -- Handle math errors


DESCRIPTION



ROUTINE DETAILS

procedure ui_math_error (err, next_handle)


ARGUMENTS
err (input: int)
error code

next_handle (input: int)
Next handler to call



UI_OPEN_DB (25Jan96) source UI_OPEN_DB (25Jan96)



NAME
ui_open_db -- Open a Fit Component Database


DESCRIPTION



ROUTINE DETAILS

procedure ui_open_db (as, file)


ARGUMENTS
as (input: pointer)
ASpec object

file[ARB] (input: char)
File name of database to open



UI_OPEN_SPEC (20Jul95) source UI_OPEN_SPEC (20Jul95)



NAME
ui_open_spec -- Open a spectrum from a file.



UI_PLOT_SPEC_WORK (25Jan96) source UI_PLOT_SPEC_WORK (25Jan96)



NAME
ui_plot_spec_work -- Plot spectra in workPlot Gterm


DESCRIPTION



ROUTINE DETAILS

procedure ui_plot_spec_work (as)


ARGUMENTS
as (input: pointer)
ASpec object



UI_PLOT_SPEC_FULL (25Jan96) source UI_PLOT_SPEC_FULL (25Jan96)



NAME
ui_plot_spec_full -- Plot spectra in fullPlot Gterm


DESCRIPTION



ROUTINE DETAILS

procedure ui_plot_spec_full (as)


ARGUMENTS
as (input: pointer)
ASpec object



UI_PLOT_AUX (25Jan96) source UI_PLOT_AUX (25Jan96)



NAME
ui_plot_aux -- Plot auxilliary information


DESCRIPTION



ROUTINE DETAILS

procedure ui_plot_aux (as)


ARGUMENTS
as (input: pointer)
ASpec object



UI_PSET (25Jan96) source UI_PSET (25Jan96)



NAME
ui_pset -- Set graphics parameters


DESCRIPTION



ROUTINE DETAILS

procedure ui_pset (pstr)


ARGUMENTS
pstr[ARB] (input: char)
parameter string



UI_REPORT (25Jan96) source UI_REPORT (25Jan96)



NAME
ui_report -- Report running information about the fitting process


DESCRIPTION



ROUTINE DETAILS

procedure ui_report (iter, chisq, flush)


ARGUMENTS
iter (input: int)
The current iteration count

chisq (input: real)
The current chisquare value

flush (input: bool)
Update the graph regardless



UI_RESIZE (25Jan96) source UI_RESIZE (25Jan96)



NAME
ui_resize -- Handle gterm resizing


DESCRIPTION



ROUTINE DETAILS

procedure ui_resize (as, gterm)


ARGUMENTS
as (input: pointer)
ASpec object

gterm[ARB] (input: char)
Gterm which is resized



UI_SAVE_DB (25Jan96) source UI_SAVE_DB (25Jan96)



NAME
ui_save_db -- Save components into a fit database


DESCRIPTION



ROUTINE DETAILS

procedure ui_save_db (as, file)


ARGUMENTS
as (input: pointer)
ASpec object

file[ARB] (input: char)
File name of database to save to



UI_SAVE_SPEC (25Jan96) source UI_SAVE_SPEC (25Jan96)



NAME
ui_save_spec -- Save a spectrum from a file


DESCRIPTION
The parameters related to saving the fit as a spectrum are as
follows:

file
The name of the table to create. This is necessary and has no
default.

(flux_col = "flux")
Name of the column containing the flux.

(wave_col = "wave")
Name of the column containing wavelength.

(sigma_col = "sigma")
Name of the column containing the flux uncertainties.

(mask_col = "mask")
Name of the column containing the data quality mask.

(scal_col = "scale")
Name of the column containing the flux scale.

(npts = 2000)
Number of points the fit should be evaluated to.


ROUTINE DETAILS

procedure ui_save_spec (as)


ARGUMENTS
as (input: pointer)
ASpec object



UI_SCROLL (25Jan96) source UI_SCROLL (25Jan96)



NAME
ui_scroll -- Set scroll bar to delineate what is in the work plot


DESCRIPTION



ROUTINE DETAILS

procedure ui_scroll ()



UI_SET_WORK (25Jan96) source UI_SET_WORK (25Jan96)



NAME
ui_set_work -- Set workplot range based on input NDC coordinates


DESCRIPTION



ROUTINE DETAILS

procedure ui_set_work (as, x1, x2, y1, y2, source, ffix)


ARGUMENTS
as (input: pointer)
ASpec object

x1 (input: double)
Range to plot, NDC relative to source

x2 (input: double)
Range to plot, NDC relative to source

y1 (input: double)
Verticle range, NDC relative to source

y2 (input: double)
Verticle range, NDC relative to source

source[ARB] (input: char)
Source which x1, x2 are relative

ffix (input: bool)
Fix Flux at current values?



UI_SET_GTERM (25Jan96) source UI_SET_GTERM (25Jan96)



NAME
ui_set_gterm -- Set the specified Gterm


DESCRIPTION



ROUTINE DETAILS

procedure ui_set_gterm (gterm)


ARGUMENTS
gterm (input: int)
Gterm ID to set to



UI_SPL_CV (25Jan96) source UI_SPL_CV (25Jan96)



NAME
ui_spl_cv -- Convert input spectrum list to UI display list


DESCRIPTION



ROUTINE DETAILS

procedure ui_spl_cv (l)


ARGUMENTS
l (input: pointer)
Internal list of spectra



UI_TO_NDC (25Jan96) source UI_TO_NDC (25Jan96)



NAME
ui_to_ndc -- Convert from User coordinates to NDC


DESCRIPTION



ROUTINE DETAILS

procedure ui_to_ndc (as, x, y, dest, nx, ny)


ARGUMENTS
as (input: pointer)
ASpec object

x (input: pointer)
NDC coordinates

y (input: pointer)
NDC coordinates

dest[ARB] (input: char)
Destination of the NDC coordinates

nx (output: double)
User coordinates

ny (output: double)
User coordinates



UI_UNDO_MSG (25Jan96) source UI_UNDO_MSG (25Jan96)



NAME
ui_undo_msg -- Inform GUI of possible undo operation


DESCRIPTION




ROUTINE DETAILS

procedure ui_undo_msg (op)


ARGUMENTS
op (input: int)
UNDO operation



UI_UNITS (25Jan96) source UI_UNITS (25Jan96)

NAME ui_units -- Get display units


DESCRIPTION



ROUTINE DETAILS

procedure ui_units (funit, dunit)


ARGUMENTS
funit (output: int)
Display flux units

dunit (output: int)
Display flux units



UI_WCS_TRANS (25Jan96) source UI_WCS_TRANS (25Jan96)



NAME
ui_wcs_trans -- Translate to/from NDC and User coordinates


DESCRIPTION



ROUTINE DETAILS

procedure ui_wcs_trans (as, dir, x, y, source)


ARGUMENTS
as (input: pointer)
ASpec object

dir[ARB] (input: char)
Direction of translation

x (input: double)
Coordinates to transform

y (input: double)
Coordinates to transform

source[ARB] (input: char)
Source/Destination of NDC coordinates



UI_WORK_SPEC (25Jan96) source UI_WORK_SPEC (25Jan96)

NAME ui_work_spec -- Determine the set of spectra in the workPlot


DESCRIPTION



ROUTINE DETAILS

pointer procedure ui_work_spec (spl)


ARGUMENTS
spl (input: pointer)
Internal spectrum list

RETURNS




UI_WORK_MARK (25Jan96) source UI_WORK_MARK (25Jan96)



NAME
ui_work_mark -- Set the marker in the fullPlot window


DESCRIPTION



ROUTINE DETAILS

procedure ui_work_mark ()


ARGUMENTS



UI_WRITE_PARAMS (25Jan96) source UI_WRITE_PARAMS (25Jan96)



NAME
ui_write_params -- Send a components parameters to the GUI


DESCRIPTION



ROUTINE DETAILS

procedure ui_write_params (c)


ARGUMENTS
c (input: pointer)
component