.. and a link to vpfit extension.

Last revision: 15 April 2003

RDGEN- for general handling of data.

This is not even vaguely supported, and is merely a collection of routines fo convenience. There is an attempt to make them  compatible with the VPFIT ones, but nothing is guaranteed. Comments here started out simply as reminders for me, and so the descriptions are pretty patchy. And if something does not work, don't be surprised ...

Routines included:

ab generate absorption list (p1:file)
ad linearize and sum spectra
bb generate a black body spectrum
cy copy data arrays
dc divide data by continuum (p1: low threshold)
ex exit program
gp generate Voigt profiles in continuum
hd dust extinction from Hb/Ha
he (or ?) give a list of commands
le set data array length
lo exit program (for IRAF devotees)
mf median filter data, result to continuum (p1: no. of pixels)
mu multiprocess files
no add noise to continuum -> new spectrum
pc close plot device
pf linelist file name input
pg plot with cursor control
pw print wavelength coefficients
qu quit program
rc recall internal data
rd read IRAF data file
sm (p1) boxcar smooth p1 pixels
sp plot the data
st store data internally
sw swap spectrum with stored one
sz data stacker
un set continuum to unity for current data
wa wavelength coefficient reset
wl set linear wavelength coefficients and length
wr write an IRAF data file
wt write data as a text (ASCII) file (p1: filename)
zr redshift correct to rest wavelength (p1: redshift, p2:n,b,[i])
zt ztable function
<  redirect input (p1: filename)


If the environment variable  RDSTART is set to a filename (like  rd_start.dat) which contains a list of commands, then these commands will be executed before the command prompt. This is useful if you always want to read in the same linelist files, or the same data. An example file is:

pf /data/cass72/rfc/qsos/
pf /data/cass72/rfc/qsos/
pf /data/cass72/rfc/qsos/
pf /data/cass72/rfc/qsos/

which will read in the four named line lists for use with cursor-controlled data plotting, and then list the available commands.

Command list:

The command list is held on a file called rd_prhelp.dat, rather than internally. There are no particular advantages, except that the program does not have to be recompiled when I remember I have not updated the help list! It goes to the file given by the environment variable RD_PRSETUP, unless you happen to be on the Institute of Astronomy, Cambridge, cluster of machines.

Data reading and moving around:

Just rd <filename> and IRAF data will be read in.

wr will write the current data to a new file. You will be prompted fo a filename if you do not give it on the command line. Note that it writes out ONLY the data, and not the error estimate or the continuum. If you want to write them as well you have to copy (cy) the data from the array in which it is held to the data array. If the file already exists, then the attempt to write it will fail.

cy gives
  Copy data array
  Enter two letter code in order (from)(to)
  continuum=c, scrunched (summed) array=s, workspace=w
  e.g. cw for cont. to work (help for more info.)

The internal arrays are workspace (w), error (e), continuum (c), scrunched (linearized summed) data (s) and its error estimate (t)
ws copies data and error to scrunched data and error
sw copies scrunched data and error to data and error
You can suppress the copying of the error arrays by appending an ''n" to the above two commands-- swn & wsn copy data but not errors
ew is error to workspace - so only one array is copied
we is workspace to error
cw is continuum to workspace
wc is workspace to continuum

A few specialist options are also available. Type "help" (after cy) to see what they are.


Looking at the data:

Plotting the data may be done using the cursor (pg) or using command line input (sp). The cursor mode operation is fairly flexible, and mimics a few of the IRAF 'splot' options. A list of one-letter commands in the graphics window is:

 Left mouse button, or "e", expands plot
 Center button, or "r" replots
 Right button, "q" or "Q" to exit
 " " print wavelength, flux at cursor position
 "." shift range up
 "," shift range down
 "=" next velocity plot to a postscript file
 "-" flag as bad data
 "a" plot whole array
 "b" mark region boundaries (B,o)
 "d" demagnify by factor 2
 "h" print line parameters for Ly-a (H)
 "j" change y of x-point nearest cursor
 "l" print line parameters (L)
 "m" wavc, me/rms, mean, rms, mean error
 "u" renormalized velocity scale
 "v" velocity scale
 "w" wavelength scale
 "x" interpolate data
 "y" max y from cursor
 "z" change continuum at point
 "*" print a *
 "<CR>" print a blank line
 "C" overplot continuum
 "G" plot data, continuum and error
 "M" mark lines
 "N" no line marking
 "P" command line prompt
 "S" suppress error, plot data only
 "Z" mark lines at reference redshift
 .. any other lower case letter for command line prompt

You can ask  (on the command line) that the continuum be plotted with a particular colour (at co co 6 will give you purple continuum), and this may contain fits to the data which are added in to the continuum using gp (see below).

One use is stacked plots of absorption lines on a common velocity scale, so you can easily see if HI, SiIV, CIV are all there together, or if low ionization lines are present, or the whole Lyman series, or low redshift Mg and Fe, or anything you care to look for. First you need to set up a file with the ion and wavelengths listed, one per line, like:

H I   1215.6701
H I   1025.7223
H I    972.5368
H I    949.7431
H I    937.8035
C III  977.020
C IV  1548.195
C IV  1550.770
N III  989.799
N V   1238.821
N V   1242.804
O VI  1031.927
O VI  1037.616
SiIII 1206.500
SiIV  1393.755
SiIV  1402.770

Then add it to an access list using pf <filename>, where <filename> is whatever you have called it. You can do this for up to 20 files, and select which you use later from a numbered list.

With this set up, you can use pg for velocity plots.
 at co   co   3
 at er   co   2
 at ti   co   5
 plot parameter? (type he for options list)
> sp
 plot parameter? (type he for options list)

 PGPLOT device? (? for list)
> xw
 Expand plot if needed: [program is now in cursor mode, so enter letter commands etc from the pgplot window]
 Cursor ("e" to mark edges, "q" when OK)
left mouse button
  Channel  19906
 right edge..
right mouse button
  Channel  24638
 Cursor ("e" to mark edges, "q" when OK) [mark off line of interest]
v    (in cursor window)
 File ID?
  .. in plot window
1    (in cursor window)

A list of possible ID's will appear on the plot. Click the left mouse button on the appropriate one and a velocity scaled plot will appear for all the lines in the list, with each marked. It should look roughly like this (click on the picture for a readable version):

You can shift the centre by moving the cursor to a new position and entering 'v' again, as often as you want. To look at a particular line in detail, just use the left cursor button to mark the left and right positions on the screen, making sure the y-position is about level with the label, and you will get it on a wavelength scale. Using 'v' and the appropriate line ID you can go back to the velocity plot if you want to. When you are fed up with this, 'a' gives you the full spectrum back -- and if you get lost (or the program does), use 'a' to restart. 'q' quits this mode, so you can read in more data or whatever.

If you are visually seaching for a system, then an easy way to do it is to set up the lines you want in a file, start with 'v' or 'u' as above, and then put the cursor near the right (for searching up in wavelength) or left (for searching down), and hit 'v' (or 'u') repeatedly until the lines you want in the pattern appear together.  In this mode the ',' and '.' keys have the same effect as 'v', so, depending on the position of the cursor, you may find yourself going in the opposite  direction in wavelength to the one you expect! The redshift of the zero velocity line in the plot is given at the top of the pgplot window.

You can send the stacked velocity plots to a postscript file, for which the best approach is to create one on the screen using 'v', when you are happy type '=', then place the cursor on the centre of the line you want and type 'v' again. The interactive screen will vanish while a hardcopy is done, and then what you sent to a file will appear on the screen. Since you might want to do this more than once, the files are named in order, etc so you don't overwrite the hardcopy each time you use '='.You should then print the plotfiles, or examine and print at leisure. If you exceed 99 in a session I have no idea what happens. If you restart the program, it starts at again, so any old one by that name will be overwritten.

You can also send other plots from the interactive plot mode to a postscript file. A way do this is to set up the display you want, then '=', and then type 'w' to redraw the plot. An alternative approach is to use the non-interactive plotting routines ('pl') and set up a list of commands which will do similar for you, and then either cut and paste into the command window or set up as a redirected input file. This is the only method of getting velocity stacked plots with data from multiple files.

If you type 'M' then wavelength (but not velocity) plots have the line positions from all ions in the input Voigt profile generating file marked, with identifications. This can be useful when unravelling messy blended systems.

'Z' takes the last reference redshift used (from 'v'), and plots the positions of all the lines in the atom.dat file at that redshift when you replot on a wavelength scale.

'N' removes the line position markers, and leaves you with an uncluttered plot.


Generating Voigt Profiles:

Voigt profiles, convolved with a Guassian instrument profile, can be inserted in the continuum array.The parameters may be entered by hand, or read from files generated by the VPFIT program.

 ion,level,col,bval,zed?    ..  or ...
 <filename> (fmt,clo,chi,zlo,zhi,blo,bhi,binc,type)

If a filename is given (so don't make SiII, or HI a filename!), then the format is needed, and you may care to restrict the ranges of the parameters for inclusion:
clo < log column density < chi,
zlo < redshift < zhi, and
blo < Doppler parameter < bhi  are included, and anything with
Doppler parameter < binc whether or not the other criteria are satisfied.
The final parameter in the list, 'type', if present, is used to include or exclude classes of ions. If you put 'l' there, only Lyman alpha will be put in for each HI found in the list. It will also suppress the Lyman limit for this or any subsequent profile generation. There are other parameters you can put in here: 'i' results in the special parameters '<>' and '__' being ignored, and 'n' results in metals from the list being omitted. You can string these together in any order, so 'nil' results in Ly-a only, with no metals and no special parameters.

Special parameters, such as <> and __, are restricted in range to the wavelength limits of the regions (marked by %% in the 26-style list) which contain them.

If you want to save the result, then you can either write out an ASCII file  (wt filename) or copy to data area and then write as IRAF .imh (cy cw followed by  wr filename).

Note added 18.12.02:
There is one incompatibility between the way things are done here and for vpfit, and that concerns the interpolation between the last Lyman series line in the list and the Lyman limit. In VPFIT this is done properly, but here it is just interpolation between the residual flux at the last line and the limit, and if other lines were put in earlier this gives the wrong result. The easy way around this is to put systems with appreciable Lyman limits in first in the input file, then they will be OK. It might be fixed sometime.


Stacking spectra:

This is a multistage process at present -- you need to fit Voigt profiles, divide by the continuum, and then stack the redshift corrected spectra. So a combination of gp, dc, and sz is used.

Stack redshift-corrected spectra:

 <filename> (fmt,elm,ion,clo,chi,zlo,zhi,blo,bhi,binc)

The filename with the line parameters must be given, and the other parameters are optional.

Default values for parameter line:
Quantity Default Comments
fmt 26 File format. Alternative is 13. For any other value 26 is used
elm H Any from atom.dat list
ion I Any from atom.dat list which is compatible with elm
clo 0.0 Minimum log column density for inclusion
chi Infinity Maximum log column density for inclusion
zlo -1000 Minimum redshift for inclusion
zhi 1E25 Maximum redshift
blo 0.0 Minimum b-value (km/s)
bhi 1E25 Maximum b-value
binc 0.0 Include ALL lines with b-value < this, even if other criteria not satisfied


Adding noise for simulated spectra:

To add roughly Gaussian noise to the continuum just type

 random number seed
 (<CR> for one based on the internal clock)
 s/n relative to peak continuum? <CR> = infinity
 (OR c1,c2 gives sigma=sqrt(cont*c1**2 + c2**2))
 work arrays replaced with synthetic spectrum

If you want to add noise to the data array just type ' no da' instead of 'no'.


Cross-correlate spectrum with wavelength table:

In an effort to help identify redshift systems a crude cross-correlation between the deviations below the continuum and a line list of the type (above) can be invoked by 'zt'. The data is replaced by the cross-correlation as a function of redshift. It was written in the hope of extracting information on metal systems hiding in the Ly-a forest, but is not particularly useful.


FLGEN- for generating a file list input for Lyman alpha/beta fits, starting with a list of
fits based on assuming Ly-a only for everything.

To use this you should have fitted everything as Ly-a and merged the resultant f26.wxxx files into a single file. The program then seeks Ly-b regions corresponding to Ly-a with a column density above some threshold (you specify), and writes out a start file (fort.13) for refitting Ly-a/Ly-b simultaneously. There is an option to fix the Ly-a's in the Ly-b regions which are unaffected (to first order) by the Ly-a's in the chosen Ly-a region. All it does is fix everything outside the wavelength range corresponding to the Ly-b projection of the Ly-a window for this.