4._LINEAR_SPECTRAL_WORLD_COORDINATE_SYSTEMS · 5._MULTISPEC_SPECTRAL_WORLD_COORDINATE_SYSTEM

5.1_LINEAR_AND_LOG_LINEAR_DISPERSION_FUNCTION · 5.2_CHEBYSHEV_POLYNOMIAL_DISPERSION_FUNCTION

5.3_LEGENDRE_POLYNOMIAL_DISPERSION_FUNCTION · 5.4_LINEAR_SPLINE_DISPERSION_FUNCTION

5.5_CUBIC_SPLINE_DISPERSION_FUNCTION · 5.6_PIXEL_ARRAY_DISPERSION_FUNCTION

5.7_SAMPLED_ARRAY_DISPERSION_FUNCTION

**The IRAF/NOAO Spectral World Coordinate Systems**

## 1. Types of Spectral Data

Spectra are stored as one, two, or three dimensional images with one axis being the dispersion axis. A pixel value is the flux over some interval of wavelength and position. The simplest example of a spectrum is a one dimensional image which has pixel values as a function of wavelength.

There are two types of higher dimensional spectral image formats. One type
has spatial axes for the other dimensions and the dispersion axis may be
along any of the image axes. Typically this type of format is used for
long slit (two dimensional) and Fabry-Perot (three dimensional) spectra.
This type of spectra is referred to as *spatial*
spectra and the
world coordinate system (WCS) format is called *ndspec*
.
The details of the world coordinate systems are discussed later.

The second type of higher dimensional spectral image consists of multiple,
independent, one dimensional spectra stored in the higher dimensions with
the first image axis being the dispersion axis; i.e. each line is a
spectrum. This format allows associating many spectra and related
parameters in a single data object. This type of spectra is referred to
as *multispec*
and the there are two coordinate system formats,
*equispec*
and *multispec*
. The *equispec*
format applies
to the common case where all spectra have the same linear dispersion
relation. The *multispec*
format applies to the general case of spectra
with differing dispersion relations or non-linear dispersion functions.
These multi-spectrum formats are important since maintaining large numbers
of spectra as individual one dimensional images is very unwieldy for the
user and inefficient for the software.

Examples of multispec spectral images are spectra extracted from a
multi-fiber or multi-aperture spectrograph or orders from an echelle
spectrum. The second axis is some arbitrary indexing of the spectra,
called *apertures*
, and the third dimension is used for
associated quantities. The IRAF **apextract**
package may produce
multiple spectra from a CCD image in successive image lines with an
optimally weighted spectrum, a simple aperture sum spectrum, a background
spectrum, and sigma spectrum as the associated quantities along the third
dimension of the image.

Many **onedspec**
package tasks which are designed to operate on
individual one dimensional spectra may operate on spatial spectra by
summing a number of neighboring spectra across the dispersion axis. This
eliminates the need to "extract" one dimensional spectra from the natural
format of this type of data in order to use tasks oriented towards the
display and analysis of one dimensional spectra. The dispersion axis is
either given in the image header by the keyword DISPAXIS or the package
*dispaxis*
parameter. The summing factors across the
dispersion are specified by the *nsum*
package parameter.
See "help onedspec.package" for information on these parmaeters.

One dimensional spectra, whether from multispec spatial spectra, have
several associated quantities which may appear in the image header as part
of the coordinate system description. The primary identification of a
spectrum is an integer aperture number. This number must be unique within
a single image. There is also an integer beam number used for various
purposes such as discriminating object, sky, and arc spectra in
multi-fiber/multi-aperture data or to identifying the order number in
echelle data. For spectra summed from spatial spectra the aperture number
is the central line, column, or band. In 3D images the aperture index
wraps around the lowest non-dispersion axis. Since most one dimensional
spectra are derived from an integration over one or more spatial axes, two
additional aperture parameters record the aperture limits. These limits
refer to the original pixel limits along the spatial axis. This
information is primarily for record keeping but in some cases it is used
for spatial interpolation during dispersion calibration. These values are
set either by the **apextract**
tasks or when summing neighboring vectors
in spatial spectra.

An important task to be aware of for manipulating spectra between image
formats is **scopy**
. This task allows selecting spectra from multispec
images and grouping them in various ways and also "extracts" apertures from
long slit and 3D spectra simply and without resort to the more general
**apextract**
package.

## 2. World Coordinate Systems

IRAF images have three types of coordinate systems. The pixel array
coordinates of an image or image section, i.e. the lines and
columns, are called the *logical*
coordinates. The logical coordinates of
individual pixels change as sections of the image are used or extracted.
Pixel coordinates are tied to the data, i.e. are fixed to features
in the image, are called *physical*
coordinates. Initially the logical
and physical coordinates are the equivalent but differ when image sections
or other tasks which modify the sampling of the pixels are applied.

The last type of coordinate system is called the *world*
coordinate
system. Like the physical coordinates, the world coordinates are tied to
the features in the image and remain unchanged when sections of the image
are used or extracted. If a world coordinate system is not defined for an
image, the physical coordinate system is considered to be the world
coordinate system. In spectral images the world coordinate system includes
dispersion coordinates such as wavelengths. In many tasks outside the
spectroscopy packages, for example the **plot**
, **tv**
and
**images**
packages, one may select the type of coordinate system to be
used. To make plots and get coordinates in dispersion units for spectra
with these tasks one selects the "world" system. The spectral tasks always
use world coordinates.

The coordinate systems are defined in the image headers using a set of
reserved keywords which are set, changed, and updated by various tasks.
Some of the keywords consist of simple single values following the FITS
convention. Others, the WAT keywords, encode long strings of information,
one for each coordinate axis and one applying to all axes, into a set of
sequential keywords. The values of these keywords must then be pasted
together to recover the string. The long strings contain multiple pieces
called WCS *attributes*
. In general the WCS keywords should be left to
IRAF tasks to modify. However, if one wants modify them directly some
tasks which may be used are **hedit**
, **hfix**
, **wcsedit**
,
**wcsreset**
, **specshift**
, **dopcor**
, and **sapertures**
. The
first two are useful for the simple keywords, the two "wcs" tasks are
useful for the linear ndspec and equispec formats, the next two are for the
common cases of shifting the coordinate zero point or applying a doppler
correction, and the last one is the one to use for the more complex
multispec format attributes.

## 3. Physical Coordinate System

The physical coordinate system is used by the spectral tasks when there is no dispersion coordinate information (such as before dispersion calibration), to map the physical dispersion axis to the logical dispersion axis, and in the multispec world coordinate system dispersion functions which are defined in terms of physical coordinates.

The transformation between logical and physical coordinates is defined by the header keywords LTVi, LTMi_j (where i and j are axis numbers) through the vector equation

l = |m| * p + v

where l is a logical coordinate vector, p is a physical coordinate vector, v is the origin translation vector specified by the LTV keywords and |m| is the scale/rotation matrix specified by the LTM keywords. For spectra rotation terms (nondiagonal matrix elements) generally do not make sense (in fact many tasks will not work if there is a rotation) so the transformations along each axis are given by the linear equation

where l is a logical coordinate vector, p is a physical coordinate vector, v is the origin translation vector specified by the LTV keywords and |m| is the scale/rotation matrix specified by the LTM keywords. For spectra a rotation term (nondiagonal matrix elements) generally does not make sense (in fact many tasks will not work if there is a rotation) so the transformations along each axis are given by the linear equation

li = LTMi_i * pi + LTVi.

If all the LTM/LTV keywords are missing they are assumed to have zero values except that the diagonal matrix terms, LTMi_i, are assumed to be 1. Note that if some of the keywords are present then a missing LTMi_i will take the value zero which generally causes an arithmetic or matrix inversion error in the IRAF tasks.

The dimensional mapping between logical and physical axes is given by the keywords WCSDIM and WAXMAP01. The WCSDIM keyword gives the dimensionality of the physical and world coordinate system. There must be coordinate information for that many axes in the header (though some may be missing and take their default values). If the WCSDIM keyword is missing it is assumed to be the same as the logical image dimensionality.

The syntax of the WAXMAP keyword are pairs of integer values,
one for each physical axis. The first number of each pair indicates which
current *logical*
axis corresponds to the original *physical*
axis
(in order) or zero if that axis is missing. When the first number is zero
the second number gives the offset to the element of the original axis
which is missing. As an example consider a three dimensional image in
which the second plane is extracted (an IRAF image section of [*,2,*]).
The keyword would then appear as WAXMAP01 = `1 0 0 1 2 0`. If this keyword
is missing the mapping is 1:1; i.e. the dimensionality and order of the
axes are the same.

The dimensional mapping is important because the dispersion axis for the nspec spatial spectra as specified by the DISPAXIS keyword or task parameter, or the axis definitions for the equispec and or multispec formats are always in terms of the original physical axes.

## 4. Linear Spectral World Coordinate Systems

When there is a linear or logarithmic relation between pixels and dispersion coordinates which is the same for all spectra the WCS header format is simple and uses the FITS convention (with the CD matrix keywords proposed by Hanisch and Wells 1992) for the logical pixel to world coordinate transformation. This format applies to one, two, and three dimensional data. The higher dimensional data may have either linear spatial axes or the equispec format where each one dimensional spectrum stored along the lines of the image has the same dispersion.

The FITS image header keywords describing the spectral world coordinates are CTYPEi, CRPIXi, CRVALi, and CDi_j where i and j are axis numbers. As with the physical coordinate transformation the nondiagonal or rotation terms are not expected in the spectral WCS and may cause problems if they are not zero. The CTYPEi keywords will have the value LINEAR to identify the type of of coordinate system. The transformation between dispersion coordinate, wi, and logical pixel coordinate, li, along axis i is given by

wi = CRVALi + CDi_i * (li - CRPIXi)

If the keywords are missing then the values are assumed to be zero except for the diagonal elements of the scale/rotation matrix, the CDi_i, which are assumed to be 1. If only some of the keywords are present then any missing CDi_i keywords will take the value 0 which will cause IRAF tasks to fail with arithmetic or matrix inversion errors. If the CTYPEi keyword is missing it is assumed to be "LINEAR".

If the pixel sampling is logarithmic in the dispersion coordinate, as required for radial velocity cross-correlations, the WCS coordinate values are logarithmic and wi (above) is the logarithm of the dispersion coordinate. The spectral tasks (though not other tasks) will recognize this case and automatically apply the anti-log. The two types of pixel sampling are identified by the value of the keyword DC-FLAG. A value of 0 defines a linear sampling of the dispersion and a value of 1 defines a logarithmic sampling of the dispersion. Thus, in all cases the spectral tasks will display and analyze the spectra in the same dispersion units regardless of the pixel sampling.

Other keywords which may be present are DISPAXIS for 2 and 3 dimensional spatial spectra, and the WCS attributes "system", "wtype", "label", and "units". The system attribute will usually have the value "world" for spatial spectra and "equispec" for equispec spectra. The wtype attribute will have the value "linear". Currently the label will be either "Pixel" or "Wavelength" and the units will be "Angstroms" for dispersion corrected spectra. In the future there will be more generality in the units for dispersion calibrated spectra.

Figure 1 shows the WCS keywords for a two dimensional long slit spectrum.
The coordinate system is defined to be a generic "world" system and the
wtype attributes and CTYPE keywords define the axes to be linear. The
other attributes define a label and unit for the second axis, which is the
dispersion axis as indicated by the DISPAXIS keyword. The LTM/LTV keywords
in this example show that a subsection of the original image has been
extracted with a factor of 2 block averaging along the dispersion axis.
The dispersion coordinates are given in terms of the *logical*
pixel
coordinates by the FITS keywords as defined previously.

Figure 1: Long Slit Spectrum

WAT0_001= 'system=world' WAT1_001= 'wtype=linear' WAT2_001= 'wtype=linear label=Wavelength units=Angstroms' WCSDIM = 2 DISPAXIS= 2 DC-FLAG = 0 CTYPE1 = 'LINEAR ' LTV1 = -10. LTM1_1 = 1. CRPIX1 = -9. CRVAL1 = 19.5743865966797 CD1_1 = 1.01503419876099 CTYPE2 = 'LINEAR ' LTV2 = -49.5 LTM2_2 = 0.5 CRPIX2 = -49. CRVAL2 = 4204.462890625 CD2_2 = 12.3337936401367

Figure 2 shows the WCS keywords for a three dimensional image where each line is an independent spectrum or associated data but where all spectra have the same linear dispersion. This type of coordinate system has the system name "equispec". The ancillary information about each aperture is found in the APNUM keywords. These give the aperture number, beam number, and extraction limits. In this example the LTM/LTV keywords have their default values; i.e. the logical and physical coordinates are the same.

Figure 2: Equispec Spectrum

WAT0_001= 'system=equispec' WAT1_001= 'wtype=linear label=Wavelength units=Angstroms' WAT2_001= 'wtype=linear' WAT3_001= 'wtype=linear' WCSDIM = 3 DC-FLAG = 0 APNUM1 = '41 3 7.37 13.48' APNUM2 = '15 1 28.04 34.15' APNUM3 = '33 2 43.20 49.32' CTYPE1 = 'LINEAR ' LTM1_1 = 1. CRPIX1 = 1. CRVAL1 = 4204.463 CD1_1 = 6.16689700000001 CTYPE2 = 'LINEAR ' LTM2_2 = 1. CD2_2 = 1. CTYPE3 = 'LINEAR ' LTM3_3 = 1. CD3_3 = 1.

## 5. Multispec Spectral World Coordinate System

The *multispec*
spectral world coordinate system applies only to one
dimensional spectra; i.e. there is no analog for the spatial type spectra.
It is used either when there are multiple 1D spectra with differing
dispersion functions in a single image or when the dispersion functions are
nonlinear.

The multispec coordinate system is always two dimensional though there may be an independent third axis. The two axes are coupled and they both have axis type "multispec". When the image is one dimensional the physical line is given by the dimensional reduction keyword WAXMAP. The second, line axis, has world coordinates of aperture number. The aperture numbers are integer values and need not be in any particular order but do need to be unique. This aspect of the WCS is not of particular user interest but applications use the inverse world to physical transformation to select a spectrum line given a specified aperture.

The dispersion functions are specified by attribute strings with the
identifier *specN*
where N is the *physical*
image line. The
attribute strings contain a series of numeric fields. The fields are
indicated symbolically as follows.

specN = ap beam dtype w1 dw nw z aplow aphigh [functions_i]

where there are zero or more functions having the following fields,

function_i = wt_i w0_i ftype_i [parameters] [coefficients]

The first nine fields in the attribute are common to all the dispersion
functions. The first field of the WCS attribute is the aperture number,
the second field is the beam number, and the third field is the dispersion
type with the same function as DC-FLAG in the *nspec*
and
*equispec*
formats. A value of -1 indicates the coordinates are not
dispersion coordinates (the spectrum is not dispersion calibrated), a value
of 0 indicates linear dispersion sampling, a value of 1 indicates
log-linear dispersion sampling, and a value of 2 indicates a nonlinear
dispersion.

The next two fields are the dispersion coordinate of the first
*physical*
pixel and the average dispersion interval per *physical*
pixel. For linear and log-linear dispersion types the dispersion
parameters are exact while for the nonlinear dispersion functions they are
approximate. The next field is the number of valid pixels, hence it is
possible to have spectra with varying lengths in the same image. In that
case the image is as big as the biggest spectrum and the number of pixels
selects the actual data in each image line. The next (seventh) field is a
doppler factor. This doppler factor is applied to all dispersion
coordinates by multiplying by 1/(1+z) (assuming wavelength dispersion
units). Thus a value of 0 is no doppler correction. The last two fields
are extraction aperture limits as discussed previously.

Following these fields are zero or more function descriptions. For linear or log-linear dispersion coordinate systems there are no function fields. For the nonlinear dispersion systems the function fields specify a weight, a zero point offset, the type of dispersion function, and the parameters and coefficients describing it. The function type codes, ftype_i, are 1 for a chebyshev polynomial, 2 for a legendre polynomial, 3 for a cubic spline, 4 for a linear spline, 5 for a pixel coordinate array, and 6 for a sampled coordinate array. The number of fields before the next function and the number of functions are determined from the parameters of the preceding function until the end of the attribute is reached.

The equation below shows how the final wavelength is computed based on the nfunc individual dispersion functions W_i(p). Note that this is completely general in that different function types may be combined. However, in practice when multiple functions are used they are generally of the same type and represent a calibration before and after the actual object observation with the weights based on the relative time difference between the calibration dispersion functions and the object observation.

w = sum from i=1 to nfunc {wt_i * (w0_i + W_i(p)) / (1 + z)}

The multispec coordinate systems define a transformation between physical pixel, p, and world coordinates, w. Generally there is an intermediate coordinate system used. The following equations define these coordinates. The first one shows the transformation between logical, l, and physical, p, coordinates based on the LTM/LTV keywords. The polynomial functions are defined in terms of a normalized coordinate, n, as shown in the second equation. The normalized coordinates run between -1 and 1 over the range of physical coordinates, pmin and pmax which are parameters of the function, upon which the coefficients were defined. The spline functions map the physical range into an index over the number of evenly divided spline pieces, npieces, which is a parameter of the function. This mapping is shown in the third and fourth equations where s is the continuous spline coordinate and j is the nearest integer less than or equal to s.

p = (l - LTV1) / LTM1_1 n = (p - pmiddle) / (prange / 2) = (p - (pmax+pmin)/2) / ((pmax-pmin) / 2) s = (p - pmin) / (pmax - pmin) * npieces j = int(s)

## 5.1 Linear and Log Linear Dispersion Function

The linear and log-linear dispersion functions are described by a
wavelength at the first *physical*
pixel and a wavelength increment per
*physical*
pixel. A doppler correction may also be applied. The
equations below show the two forms. Note that the coordinates returned are
always wavelength even though the pixel sampling and the dispersion
parameters may be log-linear.

w = (w1 + dw * (p - 1)) / (1 + z) w = 10 ** {(w1 + dw * (p - 1)) / (1 + z)}

Figure 3 shows an example from a multispec image with independent linear dispersion coordinates. This is a linearized echelle spectrum where each order (identified by the beam number) is stored as a separate image line.

Figure 3: Echelle Spectrum with Linear Dispersion Function

WAT0_001= 'system=multispec' WAT1_001= 'wtype=multispec label=Wavelength units=Angstroms' WAT2_001= 'wtype=multispec spec1 = "1 113 0 4955.44287109375 0.05... WAT2_002= '5 256 0. 23.22 31.27" spec2 = "2 112 0 4999.0810546875... WAT2_003= '58854293 256 0. 46.09 58.44" spec3 = "3 111 0 5043.505... WAT2_004= '928358078002 256 0. 69.28 77.89" WCSDIM = 2 CTYPE1 = 'MULTISPE' LTM1_1 = 1. CD1_1 = 1. CTYPE2 = 'MULTISPE' LTM2_2 = 1. CD2_2 = 1.

## 5.2 Chebyshev Polynomial Dispersion Function

The parameters for the chebyshev polynomial dispersion function are the order (number of coefficients) and the normalizing range of physical coordinates, pmin and pmax, over which the function is defined and which are used to compute n. Following the parameters are the order coefficients, ci. The equation below shows how to evaluate the function using an iterative definition where x_1 = 1, x_2 = n, and x_i = 2 * n * x_{i-1} - x_{i-2}.

The parameters for the chebyshev polynomial dispersion function are the order (number of coefficients) and the normalizing range of physical coordinates, pmin and pmax, over which the function is defined and which are used to compute n. Following the parameters are the order coefficients, c_i. The equation below shows how to evaluate the function using an iterative definition where x_1 = 1, x_2 = n, and x_i = 2 * n * x_{i-1} - x_{i-2}.

W = sum from i=1 to order {c_i * x_i}

## 5.3 Legendre Polynomial Dispersion Function

The parameters for the legendre polynomial dispersion function are the order (number of coefficients) and the normalizing range of physical coordinates, pmin and pmax, over which the function is defined and which are used to compute n. Following the parameters are the order coefficients, c_i. The equation below shows how to evaluate the function using an iterative definition where x_1 = 1, x_2 = n, and x_i = ((2i-3)*n*x_{i-1}-(i-2)*x_{i-2})/(i-1).

W = sum from i=1 to order {c_i * x_i}

Figure 4 shows an example from a multispec image with independent nonlinear
dispersion coordinates. This is again from an echelle spectrum. Note that
the IRAF **echelle**
package determines a two dimensional dispersion
function, in this case a bidimensional legendre polynomial, with the
independent variables being the order number and the extracted pixel
coordinate. To assign and store this function in the image is simply a
matter of collapsing the two dimensional dispersion function by fixing the
order number and combining all the terms with the same order.

Figure 4: Echelle Spectrum with Legendre Polynomial Function

WAT0_001= 'system=multispec' WAT1_001= 'wtype=multispec label=Wavelength units=Angstroms' WAT2_001= 'wtype=multispec spec1 = "1 113 2 4955.442888635351 0.05... WAT2_002= '83 256 0. 23.22 31.27 1. 0. 2 4 1. 256. 4963.0163112090... WAT2_003= '976664 -0.3191636898579552 -0.8169352858733255" spec2 =... WAT2_004= '9.081188912082 0.06387049476832223 256 0. 46.09 58.44 1... WAT2_005= '56. 5007.401409453303 8.555959076467951 -0.176732458267... WAT2_006= '09935064388" spec3 = "3 111 2 5043.505764869474 0.07097... WAT2_007= '256 0. 69.28 77.89 1. 0. 2 4 1. 256. 5052.586239197408 ... WAT2_008= '271 -0.03173489817897474 -7.190562320405975E-4" WCSDIM = 2 CTYPE1 = 'MULTISPE' LTM1_1 = 1. CD1_1 = 1. CTYPE2 = 'MULTISPE' LTM2_2 = 1. CD2_2 = 1.

## 5.4 Linear Spline Dispersion Function

The parameters for the linear spline dispersion function are the number of spline pieces, npieces, and the range of physical coordinates, pmin and pmax, over which the function is defined and which are used to compute the spline coordinate s. Following the parameters are the npieces+1 coefficients, c_i. The two coefficients used in a linear combination are selected based on the spline coordinate, where a and b are the fractions of the interval in the spline piece between the spline knots, a=(j+1)-s, b=s-j, and x_0=a, and x_1=b.

W = sum from i=0 to 1 {c_(i+j) * x_i}

## 5.5 Cubic Spline Dispersion Function

The parameters for the cubic spline dispersion function are the number of spline pieces, npieces, and the range of physical coordinates, pmin and pmax, over which the function is defined and which are used to compute the spline coordinate s. Following the parameters are the npieces+3 coefficients, c_i. The four coefficients used are selected based on the spline coordinate. The fractions of the interval between the integer spline knots are given by a and b, a=(j+1)-s, b=s-j, and x_0 =a sup 3, x_1 =(1+3*a*(1+a*b)), x_2 =(1+3*b*(1+a*b)), and x_3 =b**3.

The parameters for the cubic spline dispersion function are the number of spline pieces, npieces, and the range of physical coordinates, pmin and pmax, over which the function is defined and which are used to compute the spline coordinate s. Following the parameters are the npieces+3 coefficients, c_i. The four coefficients used are selected based on the spline coordinate. The fractions of the interval between the integer spline knots are given by a and b, a=(j+1)-s, b=s-j, and x_0=a**3, x_1=(1+3*a*(1+a*b)), x_2=(1+3*b*(1+a*b)), and x_3=b**3.

W = sum from i=0 to 3 {c_(i+j) * x_i}

## 5.6 Pixel Array Dispersion Function

The parameters for the pixel array dispersion function consists of just the number of coordinates ncoords. Following this are the wavelengths at integer physical pixel coordinates starting with 1. To evaluate a wavelength at some physical coordinate, not necessarily an integer, a linear interpolation is used between the nearest integer physical coordinates and the desired physical coordinate where a and b are the usual fractional intervals k=int(p), a=(k+1)-p, b=p-k, and x_0=a, and x_1=b.

W = sum from i=0 to 1 {c_(i+j) * x_i}

## 5.7 Sampled Array Dispersion Function

The parameters for the sampled array dispersion function consists of the number of coordinate pairs, ncoords, and a dummy field. Following these are the physical coordinate and wavelength pairs which are in increasing order. The nearest physical coordinates to the desired physical coordinate are located and a linear interpolation is computed between the two sample points.