splot -- plot and analyze spectra
splot images [line [band]]
The following parameters are used for error estimates in the d, k, and e key measurements. See the ERROR ESTIMATES section for a discussion of the error estimates.
sigma**2 = sigma0**2 + invgain * I
where I is the pixel value and "**2" means the square of the quantity. If either parameter is specified as INDEF or with a value less than zero then no sigma estimates are made and so no error estimates for the measured parameters are made.
The following parameters are for the interactive curve fitting function entered with the t key. This function is usually used for continuum fitting. The values of these parameters are updated during the fitting. See icfit for additional details on interactive curve fitting.
The following parameters are used to overplot standard star fluxes with the y key. See standard for more information about these parameters.
The following parameters are used for queries in response to particular keystrokes.
Splot provides an interactive facility to display and analyze spectra. See also bplot for a version of this task useful for making many plots noninteractively. Each spectrum in the image list is displayed successively. To quit the current image and go on to the next the q cursor command is used. If an image is two-dimensional, such as with multiple aperture or long slit spectra, the aperture or image column/line to be displayed is needed. If the image is three-dimensional, such as with the extra information produced by apextract , the band is needed. These parameters are queried unless specified on the command line. If given on the command line it will not be possible to change them interactively.
The plots are made on the specfied graphics device which is usually to the graphics terminal. The initial plot limits are set with the parameters xmin, xmax, ymin , and ymax . If a limit is INDEF then that limit is determined from the range of the data. The "zero" option may also be set in the options parameter to set the lower intensity limit to zero. Other options that may be set to control the initial plot are to exclude the system identification banner, and to select a histogram line type instead of connecting the pixel centers. The dispersion units used in the plot are set by the units parameter. This allows converting to units other than those in which the dispersion coordinates are defined in the spectra.
The option parameter, mentioned in the previous paragraph, is a a list of zero or more options. As previously noted, some of the options control the initial appearance of the plots. The "auto" option determines how frequently plots are redrawn. For slow terminals or via modems one might wish to minimize the redrawing. The default, however, is to redraw when changes are made. The "xydraw" parameter is specific to the x key.
After the initial graph is made an interactive cursor loop is entered. The cursor parameter may be reset to read from a file but generally the graphics device cursor is read. The cursor loop takes single keystroke commands and typed in commands begun with a colon, called colon commands. These commands are described below and a summary of the commands may be produced interactively with the ? key or a scrolling help on the status line with the / key.
Modifications to the spectra being analyzed may be saved using the i key in a new, the current, or other existing spectra. A new image is created as a new copy of the current spectrum and so if the current spectrum is part of a multiple spectrum image (including a long slit spectrum) the other spectra are copied. If other spectra in the same image are then modified and saved use the overwrite option to replace then in the new output image. If the output spectrum already exists then the overwrite flag must be set to allow modifying the data. This includes the case when the output spectrum is the same as the input spectrum. The only odd case here is when the input spectrum is one dimensional and the output spectrum is two dimensional. In this case the user is queried for the line to be written.
The other form of output, apart from that produced on the terminal, are measurements of equivalent widths, and other analysis functions. This information will be recorded in the save_file if specified.
The following keystrokes are active in addition to the normal IRAF cursor facilities (available with ":.help"):
The following keystrokes are available in the function mode. Binary operations with a constant or a second spectrum produce a query for the constant value or spectrum name.
The label must be quoted if it contains blanks. A label beginning with % (i.e. %.2f) is treated as a format for the x cursor position. The optional format is a gtext string (see help on "cursors"). The labels are not remembered between redraws.
The single profile (k) and multiple profile deblending (d) commands fit gaussian, lorentzian, and voigt line profiles with a linear background. The single profile fit, k key, is a special case of the multiple profile fitting designed to be simple to use. Two cursor positions define the region to be fit and a fixed linear continuum. The second key is used to select the type of profile to fit with g for gaussian, l for lorentzian, and v for voigt. Any other second key will default to a gaussian profile. The profile center, peak strength, and width(s) are then determined and the results are printed on the status line and in the log file. The meaning of these quantities is described later. The fit is also overplotted and may be subtracted from the spectrum subsequently with the - key.
The more complex deblending function, d key, defines the fitting region and initial linear continuum in the same way with two cursor positions. The continuum may be included in the fitting as an option. The lines to be fit are entered with the cursor near the line center (g for gaussian, l for lorentzian, v for voigt), by typing the wavelengths (t), or read from a file (f). The latter two methods are useful if the wavelengths of the lines are known accurately and if fits restricting the absolute or relative positions of the lines will be used. The t key is restricted to gaussian fits only.
The f key asks for a line list file. The format of this file has one or more columns. The columns are the wavelength, the peak value (relative to the continuum with negative values being absorption), the profile type (gaussian, lorentzian, or voigt), and the gaussian and/or lorentzian FWHM. End columns may be missing or INDEF values may be used to have values be approximated. Below are examples of the file line formats
wavelength wavelength peak wavelength peak (gaussian|lorenzian|voigt) wavelength peak gaussian gfwhm wavelength peak lorentzian lfwhm wavelength peak voigt gfwhm wavelength peak voigt gfwhm lfwhm 1234.5 <- Wavelength only 1234.5 -100 <- Wavelength and peak 1234.5 INDEF v <- Wavelength and profile type 1234.5 INDEF g 12 <- Wavelength and gaussian FWHM
where peak is the peak value, gfwhm is the gaussian FWHM, and lfwhm is the lorentzian FWHM. This format is the same as used by fitprofs and also by artdata.mk1dspec (except in the latter case the peak is normalized to a continuum of 1).
There are four queries made to define the set of parameters to be fit or constrained. The positions may be held "fixed" at their input values, allowed to shift by a "single" offset from the input values, or "all" positions may be fit independently. The widths may be constrained to a "single" value or "all" fit independently. The linear background may be included in the fit or kept fixed at that input using the cursor.
As noted above, sometimes the absolute or relative wavelengths of the lines are known a priori and this information may be entered by typing the wavelengths explicitly using the t option or read from a file using the f option during marking. In this case one should fix or fit a single shift for the position. The latter may be useful if the lines are known but there is a measureable doppler shift.
After the fit, the modeled lines are overplotted. The line center, flux, equivalent width, and full width half maxima are printed on the status line for the first line. The values for the other lines and the RMS of the fit may be examined by scrolling the status line using the +, -, and r keys. To continue enter q.
The fitting may be repeated with different options until exited with q. For each line in the blend the line center, continuum intensity at the line center, the core intensity above or below the continuum, the FWHM for the gaussian and lorentzian parts, the flux above or below the continuum, and the equivalent width are recorded in the log file. All these parameters except the continuum are based on the fitted analytic profiles. Thus, even though the fitted region may not extend into the wings of a line the equivalent width measurements include the wings in the fitted profile. For direct integration of the flux use the e key.
The fitted model may be subtracted from the data (after exiting the deblending function) using the - (minus) keystroke to delimit the region for which the subtraction is to be performed. This allows you to fit a portion of a line which may be contaminated by a blend and then subtract away the entire line to examine the remaining components.
The fitting uses an interative algorithm based on the Levenberg-Marquardt method. The iterations attempt to improve the fit by varying the parameters along the gradient of improvement in the chi square. This method requires that the initial values for the parameters be close enough that the gradient leads to the correct solution rather than an incorrect local minimum in the chi square. The initial values are determined as follows:
1. If the lines are input from a data file then those values in the file are used. Missing information is determined as below. 2. The line centers are those specified by the user either by marking with the cursor, entering the wavelenths, for read from a file. 3. The initial widths are obtained by dividing the width of the marked fitting region by the number of lines and then dividing this width by a factor depending on the profile type. 4. The initial peak intensities are the data values at the given line centers with the marked continuum subtracted.
Note that each time a new fitting option is specified the initial parameters are those from the previous fits. Thus the results do depend on the history of previous fits until the fitting is exited. Within each fit an iteration of parameters is performed as described next.
The iteration is more likely to fail if one initially attempts to fit too many parameters simultaneously. A constrained approach to the solution is obtained by iterating starting with a few parameters and then adding more parameters as the solution approaches the true chi square minimum. This is done by using the solutions from the more constrained options as the starting point for the less constrained options. In particular, the positions and a single width are fit first with fixed background. Then multiple widths and the background are added.
To conclude, here are some general comments. The most restrictive (fixed positions and single width(s)) will give odd results if the initial positions are not close to the true centers. The most general (simultaneous positions, widths, and background) can also lead to incorrect results by using unphysically different widths to make one line very narrow and another very broad in an attempt to fit very blended lines. The algorithm works well when the lines are not severely blended and the shapes of the lines are close to the profile type.
There are currently five techniques in SPLOT to measure equivalent widths and other line profile parameters. The simplest (conceptually) is by integration of the pixel values between two marked pixels. This is invoked with the e keystroke. The user marks the two edges of the line at the continuum. The measured line center, contiuum value, line flux, and equivalent width are given by:
center = sum (w(i) * (I(i)-C(i))**3/2) / sum ((I(i)-C(i))**3/2) continuum = C(midpoint) flux = sum ((I(i)-C(i)) * (w(i2) - w(i1)) / (i2 - i2) eq. width = sum (1 - I(i)/C(i))
where w(i) is the wavelength of pixel i, i1 and i2 are the nearest integer pixel limits of the integrated wavelength range, I(i) is the data value of pixel i, C(i) is the continuum at pixel (i), and the sum is over the marked range of pixels. The continuum is a linear function between the two points marked. The factor mulitplying the continuum subtracted pixel values in the flux calculation is the wavelength interval per pixel so that the flux integration is done in wavelength units. (See the discussion at the end of this section concerning flux units).
The most complex method for computing line profile parameters is performed by the profile fitting and deblending commands which compute a non-linear least-squares fit to the line(s). These are invoked with the d or k keystroke. These were described in detail previously.
The fourth and fifth methods, selected with the h key, determine the equivalent width from a gaussian profile defined by a constant continuum level "cont", a core depth "core", and the width of the line "dw" at some intermediate level "Iw".
I(w) = cont + core * exp (-0.5*((w-center)/sigma)**2) sigma = dw / 2 / sqrt (2 * ln (core/Iw)) fwhm = 2.355 * sigma flux = core * sigma * sqrt (2*pi) eq. width = abs (flux) / cont
where w is wavelength.
For ease of use with a large number of lines only one cursor position is used to mark the center of the line and one flux level. Note that both the x any y cursor positions are read simultaneously. From the x cursor position the line center and core intensity are determined. The region around the specified line position is searched for a minimum or maximum and a parabola is fit to better define the extremum.
The two methods based on the simple gaussian profile model differ in how they use the y cursor position and what part of the line is used. After typing h one selects the method and whether to use the left, right, or both sides of the line by a second keystroke. The l, r, and k keys require a continuum level of one. The y cursor position defines where the width of the line is determined. The a, b, and c keys use the y cursor position to define the continuum and the line width is determined at the point half way between the line core and the continuum. In both cases the width at the appropriate level is determined by the interception of the y level with the data using linear interpolation between pixels. The one-sided measurements use the half-width on the appropriate side and the two-sided measurments use the full-width.
The adopted gaussian line profile is drawn over the spectrum and the horizontal and vertical lines show the measured line width and the depth of the line center from the continuum. This model may also be subtracted from the spectrum using the - key.
The major advantages of these methods are that only a single cursor setting (both the x and y positions are used) is required and they are fast. The l, r, and k keys give more flexibility in adjusting the width of the gaussian line at the expense or requiring that the spectrum be normalized to a unit continuum. The a, b, and c keys allow measurements at any continuum level at the expense of only using the half flux level to determine the gaussian line width.
All these methods print and record in the log file the line center, continuum intensity at the line center, the flux, and the equivalent width. For the e key the flux is directly integrated while for the other methods the fitted gaussian is integrated. In addition, for the profile fitting methods the core intensity above or below the continuum, and the FWHMs are also printed. A zero value is record for the gaussian or lorentzian width if the value is not determined by profile fit. A brief line of data for each measurement is printed on the graphics status line. To get the full output and the output from previous measurments use the command ":show". This pages the output on the text output which may involve erasing the graphics.
The integrated fluxes for all the methods are in the same units as the intensities and the integration is done in the same units as the plotted scale. It is the user's responsiblity to keep track of the flux units. As a caution, if the data is in flux per unit frequency, say ergs/cm2/sec/hz, and the dispersion in Angstroms then the integrated flux will not be in the usual units but will be A-ergs/cm2/sec/hz. For flux in wavelength units, ergs/cm2/sec/A and the dispersion scale in Angstroms the integrated flux will be correct; i.e. ergs/cm2/sec.
Note that one can compute integrated flux in pixel units by using the $ to plot in pixels. This is appropriate if the pixel values are in data numbers or photon counts to get total data number or photons.
The deblending (d), single profile fitting (k), and profile integration and equivalent width (e) functions provide error estimates for the measured parameters. This requires a model for the pixel sigmas. Currently this model is based on a Poisson statistics model of the data. The model parameters are a constant gaussian sigma and an "inverse gain" as specified by the parameters sigma0 and invgain . These parameters are used to compute the pixel value sigma from the following formula:
sigma**2 = sigma0**2 + invgain * I
where I is the pixel value and "**2" means the square of the quantity.
If either the constant sigma or the inverse gain are specified as INDEF or with values less than zero then no noise model is applied and no error estimates are computed. Also if the number of error samples is less than 10 then no error estimates are computed. Note that for processed spectra this noise model will not generally be the same as the detector readout noise and gain. These parameters would need to be estimated in some way using the statistics of the spectrum. The use of an inverse gain rather than a direct gain was choosed to allow a value of zero for this parameters. This provides a model with constant uncertainties.
The direct profile integration error estimates are computed by error propagation assuming independent pixel sigmas. Also it is assumed that the marked linear background has no errors. The error estimates are one sigma estimates. They are given in the log output (which may also be view without exiting the program using the :show command) below the value to which they apply and in parenthesis.
The deblending and profile fit error estimates are computed by Monte-Carlo simulation. The model is fit to the data (using the sigmas) and this model is used to describe the noise-free spectrum. A number of simulations, given by the nerrsample parameter, are created in which random gaussian noise is added to the noise-free spectrum using the pixel sigmas from the noise model. The model fitting is done for each simulation and the absolute deviation of each fitted parameter to model parameter is recorded. The error estimate for the each parameter is then the absolute deviation containing 68.3% of the parameter estimates. This corresponds to one sigma if the distribution of parameter estimates is gaussian though this method does not assume this.
The Monte-Carlo technique automatically includes all effects of parameter correlations and does not depend on any approximations. However the computation of the errors does take a significant amount of time. The amount of time and the accuracy of the error estimates depend on how many simulations are done. A small number of samples (of order 10) is fast but gives crude estimates. A large number (greater than 100) is slow but gives good estimates. A compromise value of 50 is recommended for many applications.
The dispersion units capability of splot allows specifying the units with the units parameter and interactively changing the units with the ":units" command. In addition the v key allows plotting in velocity units with the zero point velocity defined by the cursor position.
The units are specified by strings having a unit type from the list below along with the possible preceding modifiers, "inverse", to select the inverse of the unit and "log" to select logarithmic units. For example "log angstroms" to plot the logarithm of wavelength in Angstroms and "inv microns" to plot inverse microns. The various identifiers may be abbreviated as words but the syntax is not sophisticated enough to recognized standard scientific abbreviations except as noted below.
angstroms - Wavelength in Angstroms nanometers - Wavelength in nanometers millimicrons - Wavelength in millimicrons microns - Wavelength in microns millimeters - Wavelength in millimeters centimeter - Wavelength in centimeters meters - Wavelength in meters hertz - Frequency in hertz (cycles per second) kilohertz - Frequency in kilohertz megahertz - Frequency in megahertz gigahertz - Frequency in gigahertz m/s - Velocity in meters per second km/s - Velocity in kilometers per second ev - Energy in electron volts kev - Energy in kilo electron volts mev - Energy in mega electron volts nm - Wavelength in nanometers mm - Wavelength in millimeters cm - Wavelength in centimeters m - Wavelength in meters Hz - Frequency in hertz (cycles per second) KHz - Frequency in kilohertz MHz - Frequency in megahertz GHz - Frequency in gigahertz wn - Wave number (inverse centimeters)
The velocity units require a trailing value and unit defining the velocity zero point. For example to plot velocity relative to a wavelength of 1 micron the unit string would be:
km/s 1 micron
Some additional examples of units strings are:
milliang megahertz inv mic log hertz m/s 3 inv mic
This task has a very large number of commands and capabilities which are interactive and graphical. Therefore it these examples are fairly superficial. The user is encouraged to simply experiment with the task. To get some help use the ? or / keys.
1. To plot a single spectrum and record any measurements in the file ngc7662:
cl> splot spectrum save_file=ngc7662
2. To force all plots to display zero as the minimum y value:
cl> splot spectrum options="auto, zero"
Note that the options auto and zero can be abbreviated to one character.
3. To successively display graphs for a set of spectra with the wavelength limits set to 3000 to 6000 angstroms:
cl> splot spec* xmin=3000 xmax=6000
4. To make batch plots create a file containing the simple cursor command
0 0 0 q
or an empty file and then execute one of the following:
cl> splot spec* graphics=stdplot cursor=curfile cl> set stdvdm=splot.mc cl> splot spec* graphics=stdvdm cursor=curfile cl> splot spec* cursor=curfile >G splot.mc
The first example sends the plots to the standard plot device specified by the environment variable "stdplot". The next example sends the plots to the standard virtual display metacode file specified by the environment variable "stdvdm". The last example redirects the standard graphics to the metacode file splot.mc. To spool the metacode file the tasks stdplot and gkimosaic may be used. For a large number of plots gkimosaic is prefered since it places many plots on one page instead of one plot per page. The other GKI tasks in the plot package may be used to examine the contents of a metacode file. A simple script call bplot is provided which has the default cursor file given above and default device of "stdplot".
5. More complex plots may be produced both interactively using the = key or the ":.snap" or ":.write" commands or by preparing a script of cursor commands.
Added colon commands for labeling.
If a wavelength scale is set with either p or u then any other spectra which are not dispersion corrected will adopt this wavelength scale.
The ( and ) keys cycle through bands if there is only one spectrum.
A new option, "flip", has been added to the options parameter to select that the spectra are plotted in decreasing wavelength.
A new options "overplot" has been added to the options parameters and colon commands to permanently set overplotting. This allows quickly overplotting many spectra.
This task will now write out the current display units in the "units_display" WCS attribute. The default task units have been changed to "" to allow picking up the "units_display" units if defined.
The deblending and gaussian fitting code now subsamples the profile by a factor of 3 and fits the data pixels to the sum of the three subsamples. This accounts for finite sampling of the data.
Error estimates are provided for the deblending (d), gaussian fitting (k), and profile integration (e) results.
Because nearly every key is used there has been some shuffling, consolidating, or elimination of keys. One needs to check the run time ? help or the help to determine the key changes.
Deblending may now use any number of components and simultaneous fitting of a linear background. A new simplified version of Gaussian fitting for a single line has been added in the k key. The old k, h, and v equivalent width commands are all part of the single h command using a second key to select a specific option. The Gaussian line model from these modes may now be subtracted from the spectrum in the same way as the Gaussian fitting. The one-sided options, in particular, are interesting in this regard as a new capability.
The arithmetic functions between two spectra are now done in wavelength with resampling to a common dispersion done automatically. The t key now provides for the full power of the ICFIT package to be used on a spectrum for continuum normalization, subtraction, or line and cosmic ray removal. The x editing key may now use the nearest pixel values rather than only the y cursor position to replace regions by straight line segments. The mode is selected by the task option parameter "xydraw".
Control over the graph window (plotting limits) is better integrated so that redrawing, zooming, shifting, and the GTOOLS window commands all work well together. The new c key resets the window to the full spectrum allowing the r redraw key to redraw the current window to clean up overplots from the Gaussian fits or spectrum editing.
The dispersion units may now be selected and changed to be from hertz to Mev and the log or inverse (for wave numbers) of units taken. As part of the units package the v key or colon commands may be used to plot in velocity relative to some origin. The $ key now easily toggles between the dispersion units (whatever they may be) and pixels coordinates.
Selection of spectra has become more complex with multiaperture and long slit spectra. New keys allow selecting apertures, lines, columns, and bands as well as quickly scrolling through the lines in multiaperture spectra. Overplotting is also more general and consistent with other tasks by using the o key to toggle the next plot to be overplotted. Overplots, including those of the Gaussian line models, are now done in a different line type.
There are new colon commands to change the dispersion axis and summing parameters for 2D image, to toggle logging, and also to put comments into the log file. All the options may also be set with colon commands.
bplot, gtools, icfit, standard, package, specplot, graph, implot, fitprofs