Initial commit

11 files changed, 3351 insertions, 0 deletions

diff --git a/noao/artdata/doc/gallist.hlp b/noao/artdata/doc/gallist.hlp
new file mode 100644
index 00000000..e1f87623
--- /dev/null
+++ b/noao/artdata/doc/gallist.hlp
@@ -0,0 +1,488 @@
+.help gallist Feb90 noao.artdata
+.ih
+TASK
+gallist -- make an artificial galaxies list
+.ih
+USAGE
+gallist gallist ngals
+.ih
+PARAMETERS
+.ls gallist
+The name of the output text file for the x and y coordinates,
+magnitudes, profile types, half-flux radii, axial ratios, and position
+angles of the artificial galaxies.  Output will be appended to this
+file if it exists.
+.le
+.ls ngals = 100
+The number of galaxies in the output galaxies list.
+.le
+.ls interactive = no
+Examine plots and change the parameters of the spatial, luminosity, and
+morphology distributions interactively.
+.le
+
+			SPATIAL DISTRIBUTION
+.ls spatial = "uniform"
+Type of spatial distribution for the galaxies.  The types are:
+.ls uniform
+The galaxies are uniformly distributed between \fIxmin\fR, \fIxmax\fR,
+\fIymin\fR, and \fIymax\fR.
+.le
+.ls hubble
+The galaxies are distributed around the center of symmetry \fIxcenter\fR and
+\fIycenter\fR according to a Hubble density law of core radius
+\fIcore_radius\fR and background density \fIbase\fR.
+.le
+.ls file
+The radial density function is contained in the text file \fIsfile\fR.
+.le
+.le
+.ls xmin = 1., xmax = 512., ymin = 1., ymax = 512.
+The range of the output coordinates in pixels.
+.le
+.ls xcenter = INDEF, ycenter = INDEF
+The coordinate of the center of symmetry for the "hubble"
+and "file" radial density functions. The default is the
+midpoint of the coordinate limits.
+.le
+.ls core_radius = 50
+The core radius of the Hubble density distribution in pixels.
+.le
+.ls base = 0.0
+The background density relative to the central density of the Hubble
+density distribution.
+.le
+.ls sseed = 2
+The initial value supplied to the random number generator used to
+generate the output x and y coordinates.
+If a value of "INDEF" is given then the clock
+time (integer seconds since 1980) is used as the seed yielding
+different random numbers for each execution.
+.le
+
+			MAGNITUDE DISTRIBUTION
+.ls luminosity = "powlaw"
+Type of luminosity distribution for the galaxies.  The types are:
+.ls uniform
+The galaxies are uniformly distributed between \fIminmag\fR and
+\fImaxmag\fR.
+.le
+.ls powlaw
+The galaxies are distributed according to a power law with coefficient
+\fIpower\fR.
+.le
+.ls schecter
+The galaxies are distributed according to a Schecter luminosity
+function with characteristic magnitude \fImstar\fR and power law exponent
+\fIalpha\fR between \fIminmag\fR and \fImaxmag\fR.
+.le
+.ls file
+The luminosity function is contained in the text file \fIlfile\fR.
+.le
+.le
+.ls minmag = -7., maxmag = 0.
+The range of output relative magnitudes.
+.le
+.ls mzero = 15.
+Magnitude zero point for Schecter luminosity function.
+.le
+.ls power = 0.6
+Coefficient for the power law magnitude distribution The default value
+of 0.6 is the Euclidean value.
+.le
+.ls alpha = -1.24
+The power law exponent of the Schecter luminosity function.
+The default value is that determined by Schecter from nearby galaxies.
+.le
+.ls mstar = -21.41
+The characteristic magnitude of the Schecter luminosity function.
+.le
+.ls lseed = 2
+The initial value supplied to the random number generator used to
+generate the output magnitudes.
+If a value of "INDEF" is given then the clock
+time (integer seconds since 1980) is used as the seed yielding
+different random numbers for each execution.
+.le
+
+			MORPHOLOGY DISTRIBUTION
+.ls egalmix = 0.4
+The fraction of the galaxies that are "ellipticals" represented
+by a de Vaucouleurs surface brightness law as opposed to "spirals"
+represented by an exponential disk surface brightness law.
+.le
+.ls ar = 0.3
+Minimum elliptical galaxy axial ratio (major/minor ratio).
+.le
+.ls eradius = 20.0
+The maximum elliptical galaxy half-flux semi-major scale radius.  This is
+the radius of an elliptical galaxy with magnitude \fIminmag\fR
+before a random factor is added.  Spiral galaxies and fainter galaxies
+are scaled from this value.
+.le
+.ls sradius = 1.0
+Ratio between half-flux scale radii of spiral and elliptical models at the
+same magnitude.  For example an elliptical galaxy with magnitude
+\fIminmag\fR will have radius \fIeradius\fR while a spiral galaxy
+of the same magnitude with have radius \fIsradius\fR * \fIeradius\fR.
+.le
+.ls absorption = 1.2
+Absorption correction for edge on spirals in magnitudes.
+.le
+.ls z = 0.05
+Minimum redshift for power law distributed galaxies.  This is the
+redshift assigned galaxies of magnitude \fIminmag\fR.  The redshifts
+are assumed proportional to the square root of the apparent luminosity;
+i.e the luminosity distance proportional to redshift.  The redshift is used
+for computing the mean apparent sizes of the galaxies
+according to (1+z)**2 / z.
+.le
+
+			USER FUNCTIONS
+.ls sfile = ""
+The name of the input text file containing the sampled spatial radial
+density
+function, one sample point per line, with the radius and relative probability
+in columns one and two respectively. The sample points need not be
+uniformly spaced or normalized.
+.le
+.ls nssample = 100
+The number of points at which the spatial density function is 
+sampled. If the spatial density function is analytic or approximated
+analytically (the "hubble" option) the function is sampled
+directly. If the function is read from a file  (the "file" option) an
+initial smoothing step is performed before sampling.
+.le
+.ls sorder = 10
+The order of the spline fits used to evaluate the integrated spatial
+density function.
+.le
+.ls lfile = ""
+The name of the input text file containing the sampled luminosity
+function, one sample point per line, with the magnitude and relative
+probability in columns one and two respectively. The sample points need
+not be uniformly spaced or normalized.
+.le
+.ls nlsample = 100
+The number of points at which the luminosity function is 
+sampled. If the luminosity function is analytic or approximated
+analytically (the "uniform", "powlaw" and "schecter" options) the
+function is sampled directly.  If it is read from a file
+(the "file" option) an initial smoothing step is performed before sampling.
+.le
+.ls lorder = 10
+The order of the spline fits used to evaluate the integrated
+luminosity function.
+.le
+
+			INTERACTIVE PARAMETERS
+.ls rbinsize = 10.
+The bin size in pixels of the plotted histogram of the radial density
+distribution.
+.le
+.ls mbinsize = 0.5
+The bin size in magnitudes of the plotted histogram of the luminosity function.
+.le
+.ls dbinsize = 0.5
+The bin size in pixels of the plotted histogram of the half-power semi-major
+axis distribution.
+.le
+.ls ebinsize = 0.1
+The bin size of the plotted histogram of the axial ratio distribution.
+.le
+.ls pbinsize = 20.
+The bin size in degrees of the plotted histogram of the position angle
+distribution.
+.le
+.ls graphics = stdgraph
+The default graphics device.
+.le
+.ls cursor = ""
+The graphics cursor.
+.le
+.ih
+DESCRIPTION
+\fBGallist\fR generates a list of x and y coordinates, magnitudes,
+morphological types, half-power radii, axial ratios, and position
+angles for a sample of \fIngals\fR galaxies based on a user selected
+spatial density function \fIspatial\fR  and luminosity function
+\fIluminosity\fR and writes (appends) the results to the text file
+\fIgallist\fR. If the \fIinteractive\fR parameter is "yes" the user can
+interactively examine plots of the spatial density function, the
+radial density function,  the luminosity function, radii, axial ratios,
+and position angle distributions and alter the parameters of the task
+until a satisfactory artificial field is generated.
+
+The spatial density function generates x and y values around a center
+of symmetry defined by \fIxcenter\fR and \fIycenter\fR within the x and
+y limits \fIxmin\fR, \fIxmax\fR, \fIymin\fR and \fIymax\fR according to
+the spatial density function specified by \fIspatial\fR.  The three
+supported spatial density functions are listed below where R is the
+radial distance in pixels, P is the relative spatial density, C is a
+constant, and f is the best fitting cubic spline function to the spatial
+density function R(user), P(user) supplied by the user in the text file
+\fIsfile\fR.
+
+.nf
+  uniform:  P = C
+  hubble:   P = 1.0 / (1 + R / core_radius) ** 2 + base
+  file:     P = f (R(user), P(user))
+.fi
+
+The Hubble and user spatial density functions are sampled at
+\fInssample\fR equally spaced points, and integrated to give the
+spatial density probability function at each sampled point. The
+integrated probability function is normalized and approximated by a
+cubic spline of order \fIsorder\fR.  The x and y coordinates are
+computed by randomly sampling the integrated probability function until
+\fIngals\fR galaxies which satisfy the x and y coordinate limits
+\fIxmin\fR, \fIxmax\fR, \fIymin\fR and \fIymax\fR are generated.
+
+The luminosity function generates relative magnitude values between
+\fIminmag\fR and \fImaxmag\fR (before absorption effects are added)
+according to the luminosity function specified by \fIluminosity\fR.
+The four supported luminosity functions are listed below where M is the
+magnitude, P is the relative luminosity function, C is a constant and f
+is the best fitting cubic spline function to the luminosity function
+M(user), P(user) supplied by the user in the text file \fIlfile\fR.
+
+.nf
+  uniform:   P = C
+  powlaw:    P = C * 10. ** (power * M)
+  schecter:  P = C * 10. ** (alpha * dM) * exp (-10. ** dM)
+  file:      P = f (M(user), P(user))
+
+  where      dM = 0.4 * (mstar - M + mzero)
+.fi
+
+The uniform distribution is not very physical but may be useful for
+testing.  The power law distribution is that expected for a homogeneous
+and isotropic distribution of galaxies.  The default value of 0.6 is
+that which can be calculated simply from Euclidean geometry.  Observations
+of faint galaxies generally show a smaller value.  The Schecter
+function provides a good approximation to a galaxy cluster when
+used in conjunction with the Hubble spatial distribution (though there
+is no mass segregation applied).  The "best fit" values for the
+parameters \fImstar\fR and \fIalpha\fR are taken from the paper by
+Schecter (Ap.J 203, 297, 1976).  The \fImzero\fR parameter is used
+to convert to absolute magnitudes.  Note that it is equivalent to
+set \fImzero\fR to zero and adjust the characteristic magnitude
+to the same relative magnitude scale or to use absolute magnitudes
+directly.
+
+The Schecter and user file distributions are sampled at \fInlsample\fR
+equally spaced points, and integrated to give the luminosity
+probability function at each sampled point. The probability function is
+normalized and approximated by a cubic spline of order \fIlorder\fR.
+The magnitudes are computed by randomly sampling the integrated
+probability function until \fIngals\fR objects which satisfy the
+magnitude limits \fIminmag\fR and \fImaxmag\fR are generated.
+
+The artificial galaxies have one of two morphological types,
+"ellipticals" with a de Vaucouleurs surface brightness law and
+"spirals" with an exponential surface brightness law. The fraction
+of elliptical galaxies is set by the parameter \fIegalmix\fR.  The
+position angles of the major axis are distributed uniformly between 0.0
+and 360.0 degrees.  The axial ratio (major to minor) of the elliptical
+models is allowed to range uniformly between 1 and \fIar\fR
+(that is E0 - E7).
+
+The spiral models have inclinations, i, ranging uniformly between 0 and
+90 degrees.  The axial ratio is then given by
+
+	a/b = sqrt (sin(i)**2 * .99 + .01)
+
+which is taken from Holmberg in Galaxies and the Universe (which
+references the work of Hubble).  Note the axial ratio is limited to
+0.1 by this formula.  An internal absorption correction is then
+made based on the inclination using the relation
+
+	dM = A * (min (10, cosecant (i)) - 1) / 9
+
+where is the absorption of an edge on galaxy relative to face on and
+the cosecant is limited to 10.  Note that this correction changes
+allows galaxies with magnitudes less than \fImaxmag\fR and alters
+the luminosity function somewhat.  Or in other words, the luminosity
+function is based on absorption corrected magnitudes.
+
+The sizes of the galaxy images are scaled from the maximum half-flux
+radius of an elliptical galaxy given by the parameter \fIeradius\fR.
+This is the radius given to an elliptical galaxy of magnitude
+\fIminmag\fR (prior to adding a random factor described below).  The
+ratio between the half-flux radii of the exponential disk and de
+Vaucouleurs models at a given total magnitude is set by the parameter
+\fIsradius\fR (note this is a fraction of \fIeradius\fR and not an
+actual radius).  This allows adjusting the relative surface brightness
+of elliptical and spiral models.
+
+The distribution of sizes is based on the apparent
+magnitude of the galaxies.  For the power law magnitude distribution
+the cosmological redshift factor for angular diameters is used.  The
+redshift/magnitude relation is assumed to be such that the redshift is
+proportional to the luminosity distance (the square root of the
+apparent luminosity).  Thus,
+
+
+.nf
+                Z = z * 10. ** (0.2 * (M - minmag))
+                Zfactor = ((1+Z)**2 / Z) / ((1+z)**2 / z)
+  ellipticals:  r = eradisus * Zfactor
+  spirals:      r = sradius * eradius * Zfactor
+.fi
+
+where z is the reference redshift at the minimum magnitude, and Z is the
+redshift at magnitude M.  For very small z the size varies as the
+luminosity distance but at larger z the images appear more extended with
+lower surface brightness.  For very deep simulations a pure luminosity
+distance relation gives faint galaxies which are too small and compact
+compared to actual observations.
+
+For the other magnitude distributions, the Schecter cluster function
+in particular where all galaxies are at the same distance, the scale radius
+obeys the following relation.
+
+.nf
+  ellipticals:  r = eradius * 10. ** ((minmag - M) / 6)
+  spirals:      r = sradius * eradius * 10. ** ((minmag - M) / 6)
+.fi
+
+This relation gives the size decreasing slightly less rapidly than that
+giving a constant surface brightness.  This relation is taken from
+Holmberg (Galaxies and the Universe).
+
+A uniform random factor of 50% is added to the sizes computed for
+the power law magnitude distribution and 20% for the other distributions.
+
+The interactive spatial plot shows the positions of the galaxies, the
+galaxy type (circles are de Vaucouleurs profiles and other types are
+diamonds), and rough size.
+.ih
+CURSORS
+The following interactive keystroke commands are available from within the
+GALLIST task.
+
+.nf
+	Gallist Keystroke Commands
+
+?	Print options
+f	Fit one or more of following 
+	    Spatial density function (SDF)
+            Luminosity  function (LF)
+	    Distribution of morphological type
+	    Diameter distribution
+	    Roundness distribution
+	    Position angle distribution 
+x	Plot the x-y spatial density function
+r	Plot the histogram of the radial density function
+m	Plot the histogram of the luminosity function
+d	Plot the histogram of the diameter values
+e	Plot the histogram of the roundness values 
+p	Plot the histogram of the position angle values
+:	Colon escape commands (see below)
+q	Exit program
+.fi
+
+The following parameters can be shown or set from within the GALLIST task.
+
+.nf
+		Gallist Colon Commands
+
+:show			Show gallist parameters
+:ngal       [value]	Number of galaxies
+
+:spatial    [string]	Spatial density function (SDF) (uniform|hubble|file) 
+:xmin       [value]	Minimum X value
+:xmax       [value]	Maximum X value
+:ymin       [value]	Minimum Y value
+:ymax       [value]	Maximum Y value
+:xcenter    [value]	X center for SDF
+:ycenter    [value]	Y center for SDF
+:core       [value]	Core radius for Hubble density function
+:base       [value]	Background density for Hubble density function
+
+:luminosity [string]	Luminosity function (LF)
+			(uniform|powlaw|schecter|file)
+:minmag     [value]	Minimum magnitude
+:maxmag     [value]	Maximum magnitude
+:mzero      [value]	Magnitude zero-point of schecter LF
+:power      [value]     Power law coefficient for powlaw LF
+:alpha      [value]	Schecter parameter
+:mstar      [value]	Characteristic mag for Schecter LF
+
+:egalmix    [value]	Elliptical/Spiral galaxy ratio
+:ar         [value]     Minimum elliptical galaxy axial ratio
+:eradius    [value]     Maximum elliptical half flux radius
+:sradius    [value]     Spiral/elliptical radius at same magnitude
+:z          [value]     Minimum redshift
+:absorption [value]     Absorption correction for spirals
+
+:lfile      [string]    Name of the LF file
+:sfile	    [string]    Name of the SDF file
+:nlsample   [value]	Number of LF sample points 
+:lorder	    [value]	Order of spline approximation to the integrated LF
+:nssample   [value]	Number of SDF sample points
+:sorder	    [value]	Order of spline approximation to the integrated SDF
+
+:rbinsize   [value]	Resolution of radial SDF histogram in pixels
+:mbinsize   [value]	Resolution of magnitude histogram in magnitudes
+:dbinsize   [value]	Resolution of diameter histogram in pixels
+:ebinsize   [value]	Resolution of roundness histogram in pixels
+:pbinsize   [value]     Resolution of position angle histogram in degrees
+.fi
+.ih
+EXAMPLES
+1. Create a galaxy cluster with a power law distribution of field galaxies
+and stars as background/foreground.
+
+.nf
+    ar> gallist galaxies.dat 100 spatial=hubble lum=schecter egal=.8
+    ar> gallist galaxies.dat 500
+    ar> starlist galaxies.dat 100
+    ar> mkobjects galaxies obj=galaxies.dat gain=3 rdnoise=10 poisson+
+.fi
+
+Note that the objects are appended to the same file.  Actually making
+the image with \fBmkobjects\fR takes about 5 minutes (2.5 min cpu) on a
+SPARCstation 1.
+
+2. Examine the distributions for a uniform spatial distribution
+and power law magnitude distribution using 1000 galaxies without
+creating a data file.
+
+.nf
+    ar> gallist dev$null 1000 inter+
+	    ... an x-y plot will appear on the screen
+	    ... type r to examine the radial density function
+	    ... type m to examine the luminosity function
+	    ... type d to examine the half-flux radii distribution
+	    ... type e to examine the axial ratio distribution
+	    ... type = to make a copy of any of the plots
+	    ... type q to quit
+.fi
+.ih
+REVISIONS
+.ls GALLIST V2.11+
+The random number seeds can be set from the clock time by using the value
+"INDEF" to yield different random numbers for each execution.
+.le
+.ls GALLIST V2.11
+The default value for the minimum elliptical galaxy axial ratio was
+change to 0.3 to cover the range E0-E7 uniformly.
+.le
+.ih
+BUGS
+This is a first version and is not intended to produce a full model
+of galaxy fields.  Some of the relations used are empirical and
+simple minded with the aim being to produce reasonably realistic images.
+
+The spline approximation to the spatial density and luminosity
+probability functions can cause wiggles in the output spatial density
+and luminosity functions. Users can examine the results interactively
+and experiment with the spline order and number of sample points if
+they are not satisfied with the results of GALLIST. The default setup
+of 10 sample points per spline piece is generally satisfactory for the
+spatial density and luminosity functions supplied here.
+.ih
+SEE ALSO
+starlist mkobjects
+.endhelp
diff --git a/noao/artdata/doc/mk1dspec.hlp b/noao/artdata/doc/mk1dspec.hlp
new file mode 100644
index 00000000..b96684a8
--- /dev/null
+++ b/noao/artdata/doc/mk1dspec.hlp
@@ -0,0 +1,355 @@
+.help mk1dspec Jul95 noao.artdata
+.ih
+NAME
+mk1dspec -- Make/add artificial 1D spectra
+.ih
+USAGE
+mk1dspec input
+.ih
+PARAMETERS
+.ls input
+Spectra to create or modify.
+.le
+.ls output = ""
+Output spectra when modifying input spectra.  If no output spectra are
+given then existing spectra in the input list are modified directly.
+If an output list is given then it must match in number the input list.
+.le
+.ls ap = 1
+Image line to be created or modified in images of dimension greater than 1.
+.le
+.ls rv = 0.
+Radial velocity (km/s) or redshift, as selected by the parameter \fIz\fR,
+applied to line positions and continuum.  Velocities are converted to
+redshift using the relativistic relation 1+z = sqrt ((1+rv/c)/(1-rv/c)).
+Note the shift is not a shift in the dispersion parameters but in the
+underlying artificial spectrum.
+.le
+.ls z = no
+Is the velocity parameter a radial velocity or a redshift?
+.le
+
+WHEN CREATING NEW SPECTRA
+.ls title = ""
+Image title to be given to the spectra.  Maximum of 79 characters.
+.le
+.ls ncols = 512
+Number of columns.
+.le
+.ls naps = 1
+Number of lines or apertures.
+.le
+.ls header = "artdata$stdheader.dat"
+Image or header keyword data file.  If an image is given then the image header
+is copied.  If a file is given then the FITS format cards are copied.
+This only applies to new images.   The data file consists of lines
+in FITS format with leading whitespace ignored.  A FITS card must begin
+with an uppercase/numeric keyword.  Lines not beginning with a FITS
+keyword such as comments or lower case are ignored.  The user keyword
+output of \fBimheader\fR is an acceptable data file.  See \fBmkheader\fR
+for further information.
+.le
+.ls wstart = 4000., wend = 8000.
+Starting and ending wavelengths in Angstroms.  The dispersion is
+determined by these values and the number of columns.
+.le
+
+CONTINUUM PARAMETERS
+.ls continuum = 1000., slope = 0.
+Continuum of the starting wavelength at rest and the slope of the continuum.
+.le
+.ls temperature = 5700.
+Blackbody continuum temperature in Kelvin.  A value of 0 is used if
+no blackbody continuum is desired.  The intensity level is set by
+scaling to the continuum level of the starting wavelength at rest.
+.le
+.ls fnu = no
+Compute the continuum as flux per unit frequency (F-nu) if yes or flux per
+unit wavelength (F-lambda) if no.
+.le
+
+
+LINE PARAMETERS
+.ls lines = ""
+List of spectral line files.  Spectral line files contain lines of rest
+wavelength, peak, profile type, and widths (see the DESCRIPTION
+section).  The latter parameters may be missing or INDEF in which case they
+default to the task \fIpeak\fR, \fIprofile\fR, \fIgfwhm\fR, and \fIlfwhm\fR
+parameters (note that the \fIpeak\fR parameter is not a constant but the
+random number scaling).  If no file or a new (nonexistent) file is
+specified then a number of random lines given by the parameter \fInlines\fR
+is generated.  If a new file name is specified then the lines generated are
+recorded in the file.  If the list of spectral line files is shorter than
+the list of input spectra, the last spectral line list file is reused.
+.le
+.ls nlines = 0
+If no spectral line file or a new file is specified then the task will
+generate this number of random spectral lines.  The rest wavelengths are
+uniformly random within the limits of the spectrum, the peaks are uniformly
+random between zero and the value of the \fIpeak\fR parameter, the profile
+type is given by \fIprofile\fR, and the widths are fixed at the values of
+the \fIgfhwm\fR ad \fIlfwhm\fR parameters.  If a redshift is applied the
+rest wavelengths are shifted and repeated periodically.
+.le
+.ls profile = "gaussian" (gaussian|lorentzian|voigt)
+The default profile type for random lines or when not specified in the
+spectral line file.  The profile types are:
+
+.nf
+      gaussian - Gaussian profile
+    lorentzian - Lorentzian profile
+         voigt - Voigt profile
+.fi
+.le
+.ls peak = -0.5
+The maximum spectral line peak value when generating random lines or
+when the peak is missing from the spectral line file.
+This value is relative to the continuum unless the continuum is zero.
+Negative values are absorption lines and positive values are emission lines.
+.le
+.ls gfwhm = 20., lfwhm = 20.
+The default gaussian and lorentzian full widths at half maximum (FWHM), in
+Angstroms, used when generating random lines or when the widths are missing
+from the spectral line file.
+.le
+.ls seed = 1
+Random number seed.  If a value of "INDEF" is given then the clock
+time (integer seconds since 1980) is used as the seed yielding
+different random numbers for each execution.
+.le
+
+.ls comments = yes
+Include comments recording task parameters in the image header?
+.le
+
+PACKAGE PARAMETERS
+.ls nxsub = 10
+Number of pixel subsamples used in computing the gaussian spectral line
+profiles.
+.le
+.ls dynrange = 100000.
+The gaussian line profiles extend to infinity so a dynamic range, the ratio
+of the peak intensity to the cutoff intensity, is imposed to cutoff
+the profiles.
+.le
+.ih
+DESCRIPTION
+This task creates or modifies one dimensional spectra.  with a combination
+of blackbody and linear sloped continuum and emission and absorption
+spectral lines.  The spectral lines may be gaussian, lorentzian, or voigt
+profiles.  A velocity shift may be applied to the underlying artificial
+spectrum which is shifted into the specified observed wavelength region.
+No noise is included but may be added with the task \fBmknoise\fR.  New
+spectra are created with the specified number of pixels, wavelength range,
+and real datatype.  When \fInlines\fR is greater than 1 then an image with
+the specified number of lines is created though only the line given by the
+\fIap\fR is will have a spectrum.  Existing spectra may be modified in
+place or new spectra output.  Spectra are modified by adding the continuum
+and lines defined by the parameters.
+
+For new images a set of header keywords may be added by specifying an image
+or data file with the \fIheader\fR parameter (see also \fBmkheader\fR).  If
+a data file is specified lines beginning with FITS keywords are entered in
+the image header.  Leading whitespace is ignored and any lines beginning
+with words having lowercase and nonvalid FITS keyword characters are
+ignored.  In addition to this optional header, parameters for the
+wavelength coordinates are defined.  Finally, comments may be added to the
+image header recording the task parameters and any information from the
+line file which are not line definitions.
+
+Initially all spectra are created without a dispersion function; i.e.
+pixel coordinates.  For multiple spectra in an image this task must be
+executed for each image line to set the dispersion function and add data.
+When an image line is selected if it has a defined dispersion function that
+is used otherwise the task wavelength parameters are used.
+
+A continuum is defined by the value at the starting wavelength at rest, a
+slope, and a blackbody function of a given temperature.  The blackbody
+function is scaled to have the specified continuum value at the starting
+wavelength at rest.  The blackbody flux units are per unit wavelength
+(F-lambda).  A zero continuum value or a zero temperature will not produce a
+blackbody continuum.
+
+Spectral lines are modeled by gaussian, lorentzian, or voigt profiles of
+specified wavelength, peak, and widths.  The lines are defined in a
+spectral line file or generated randomly.  A spectral line file consists of
+text lines giving rest wavelength, peak, profile type, gaussian full width
+at half maximum and/or lorentzian full width at half maximum.  Only the
+wavelength is required and subsequent fields may be missing or given as
+INDEF.  The following table shows the possible formats where wavelength,
+peak,  gfwhm, and lfwhm are values of wavelength, peak, gaussian FWHM, and
+lorentzian FWHM.  The profile types are as shown though they may be
+abbreviated to one character.
+
+.nf
+	wavelength
+	wavelength peak
+	wavelength peak gaussian
+	wavelength peak gaussian gfwhm
+	wavelength peak gaussian gfwhm
+	wavelength peak lorentzian
+	wavelength peak lorentzian lfwhm
+	wavelength peak lorentzian lfwhm
+	wavelength peak voigt
+	wavelength peak voigt gfwhm
+	wavelength peak voigt gfwhm lfwhm
+	wavelength peak voigt gfwhm lfwhm
+.fi
+
+When a field is missing or INDEF the values given by the parameters
+\fIpeak\fR, \fIprofile\fR, \fIgfwhm\fR, and \fIlfwhm\fR are used.  If a
+peak value is missing, random values between zero and the \fIpeak\fR value
+are generated.  Note that to get random line intensities with some
+specified profile type and widths the value INDEF would be used for
+the peak field.
+
+If no spectral line file is specified or a new (nonexistent) file is named
+then the number of random lines given by the parameter \fInlines\fR is
+generated.  The rest wavelengths are uniformly random within the wavelength
+range of the spectrum and extend periodically outside this range in the
+case of an applied velocity shift, the peaks are uniformly random between
+zero and the \fIpeak\fR parameter, and the profile type and widths are
+given by the \fIprofile\fR, \fIgfwhm\fR, and \fIlfwhm\fR parameters.  If a
+new file is named then the parameters of the generated lines will be
+output.
+
+The peak values are taken relative to a positive continuum.  In other
+words the generated line profile is multiplied by the continuum (with a
+minimum of zero for fully saturated absorption lines).  If the
+continuum is less than or equal to zero, as in the case of an
+artificial arc spectrum or pure emission line spectrum, then the peak
+values are absolute intensities.  Positive peak values produce emission
+lines and negative values produce absorption lines.  Odd results will
+occur if the continuum has both positive and zero or negative values.
+
+The underlying rest spectrum may be shifted.  This is used primarily for
+testing radial velocity measuring algorithms and is not intended as a
+complete model of redshift effects.  The starting and ending wavelengths
+are not changed by redshifting; these are the instrumental observed
+wavelengths.  Input line wavelengths are specified at rest and then
+shifted into or out of the final spectrum.  To be realistic the line
+list should include wavelengths over a great enough range to cover
+all desired redshifts.  The peaks and widths are also appropriately
+modified by a redshift.  As an example, if the redshift is 1 the
+lines will appear broader by a factor of 2 and the peaks will be down
+by a factor of 2 in order to maintain the same flux.
+
+The random line generation is difficult in that one wants to have the
+same set of lines (for a given seed) observed at different redshifts.
+What is done is that the specified number of random lines is generated
+within the observed wavelength interval taken at rest.  This set is
+then repeated periodical over all wavelengths.  A redshift will then
+shift these rest lines in to or out of the observed spectrum.  If the
+lines are output, they are given at rest.  \fBNote that this
+periodicity may be important in interpreting cross correlation redshift
+tests for large shifts between template and object spectra.\fR
+
+The definitions of the continuum are also affected by a redshift.
+The reference point for the continuum level, slope, and blackbody
+continuum is the starting wavelength taken at rest.  Shifts will then
+modify the continuum level at the first pixel appropriately.  In
+particular a large redshift will shift the blackbody in such a way that
+the flux is still given by the \fIcontinuum\fR parameter at the starting
+wavelength at rest.
+.ih
+EXAMPLES
+1. Create a simple blackbody continuum between the default wavelengths.
+
+.nf
+	cl> mk1dspec bb title=Blackbody
+.fi
+
+2. Create a random absorption spectrum on a blackbody continuum without
+saving the line list.
+
+.nf
+	cl> mk1dspec bbab title=Absorption nlines=100
+.fi
+
+3. Create a random absorption spectrum with noise and cosmic rays.
+
+.nf
+	cl> mk1dspec bbab title=Absorption nlines=100
+	cl> mknoise bbab rdnoise=10 poisson+ ncos=5 energy=1000
+.fi
+
+4. Create a random emission spectrum on a blackbody continuum and save
+the line list.
+
+.nf
+	cl> mk1dspec bbem title=Emission nl=30 peak=0.6 lines=bbem.dat
+.fi
+
+5. Create an artificial random arc line spectrum.
+
+.nf
+	cl> mk1dspec arc title="Arc lines" cont=0 peak=500 nl=30
+.fi
+
+6. Create a test spectrum with a line list.
+
+.nf
+	cl> type linelist
+	4100 -.1 g 20
+	4200 -2. g 20
+	4300 -.3 g 20
+	5100 -.9 g 2
+	5200 -.9 g 4
+	5300 -.9 g 8
+	6700 .9 g 8
+	6800 .9 g 2
+	6900 .9 g 4
+	7700 .3 g 20
+	7800 .2 g 20
+	7900 .1 g 20
+	cl> mk1dspec testspec title=Test cont=500 temp=0 lines=linelist
+.fi
+
+7. Add absorption lines to a spectrum.
+
+.nf
+	cl> mk1dspec bb out=artspec cont=0 lines=STDIN
+	4300 -60
+	5000 -200
+	[EOF]
+.fi
+
+Normally the input spectrum would be a real spectrum.
+
+8. Make two spectra taken from the same set of random lines but differing
+in redshift.
+
+.nf
+	cl> mk1dspec restspec nl=30
+	cl> mk1dspec redspec rv=3000 nl=30
+	cl> mk1dspec bluespec rv=-.01 z+ nl=30
+.fi
+
+9. Make a multispec image with 5 apertures and a range of redshifts.
+
+.nf
+	cl> mk1dspec spec.ms ap=1 nl=30 rv=0 naps=5
+	cl> mk1dspec spec.ms ap=2 nl=30 rv=1000
+	cl> mk1dspec spec.ms ap=3 nl=30 rv=2000
+	cl> mk1dspec spec.ms ap=4 nl=30 rv=3000
+	cl> mk1dspec spec.ms ap=5 nl=30 rv=4000
+.fi
+.ih
+REVISIONS
+.ls MK1DSPEC V2.11+
+The random number seed can be set from the clock time by using the value
+"INDEF" to yield different random numbers for each execution.
+.le
+.ls MK1DSPEC V2.11
+Lorentzian and Voigt profiles were added and the parameters and input
+line list format were changed.  The widths are now FWHM instead of
+gaussian sigmas.
+.le
+.ls MK1DSPEC V2.10.3
+The format parameter was eliminated and the task updated to produce the
+current coordinate system format.
+.le
+.ih
+SEE ALSO
+mknoise, mk2dspec, mkheader, onedspec.sinterp
+.endhelp
diff --git a/noao/artdata/doc/mk2dspec.hlp b/noao/artdata/doc/mk2dspec.hlp
new file mode 100644
index 00000000..44336400
--- /dev/null
+++ b/noao/artdata/doc/mk2dspec.hlp
@@ -0,0 +1,207 @@
+.help mk2dspec Aug90 noao.artdata
+.ih
+NAME
+mk2dspec -- Make/add 2D spectra using 1D spectra templates
+.ih
+USAGE
+mk2dspec input
+.ih
+PARAMETERS
+.ls input
+Spectra to create or modify.
+.le
+.ls output = ""
+Output spectra when modifying input spectra.  If no output spectra are
+given then existing spectra in the input list are modified directly.
+If an output list is given then it must match in number the input list.
+.le
+.ls models = ""
+List of model parameter files.  If the list of model files is shorter than the
+list of input images then the last model file is reused.  The model
+parameter files contain lines giving one dimensional spectrum template
+name, intensity scale, type of cross dispersion profile, profile
+width in the center line, change of width per line, profile position
+in the center line, and change of position per line (see the DESCRIPTION
+section).
+.le
+.ls comments = yes
+Include comments recording task parameters in the image header?
+.le
+
+WHEN CREATING NEW SPECTRA
+.ls title = ""
+Image title to be given to the spectra.  Maximum of 79 characters.
+.le
+.ls ncols = 100, nlines = 512
+Number of columns and lines.
+.le
+.ls header = "artdata$stdheader.dat"
+Image or header keyword data file.  If an image is given then the image header
+is copied.  If a file is given then the FITS format cards are copied.
+This only applies to new images.   The data file consists of lines
+in FITS format with leading whitespace ignored.  A FITS card must begin
+with an uppercase/numeric keyword.  Lines not beginning with a FITS
+keyword such as comments or lower case are ignored.  The user keyword
+output of \fBimheader\fR is an acceptable data file.  See \fBmkheader\fR
+for further information.
+.le
+.ih
+DESCRIPTION
+This task creates or modifies two dimensional spectra by taking one
+dimensional spectra, convolving them with a spatial profile across the
+dispersion, and adding them into two dimensional images.  The one
+dimensional spectra may be real data or artificial data created with
+the task \fBmk1dspec\fR.  No noise is included but may be added with
+the task \fBmknoise\fR.  The spatial profile is fully subsampled and
+may vary in width and position along the dispersion axis.  The spatial
+axis is along the first dimension and the dispersion is along the
+second dimension.
+
+For new images a set of header keywords may be added by specifying an
+image or data file with the \fIheader\fR parameter (see also \fBmkheader\fR).
+If a data file is specified lines beginning with FITS keywords are
+entered in the image header.  Leading whitespace is ignored and any
+lines beginning with words having lowercase and nonvalid FITS keyword
+characters are ignored.  In addition, comments may be added to
+the image header recording the model file name and the contents of the
+model file.
+
+The spatial profile models are specified in one or more model parameter
+files.  These files contain lines giving a one dimensional spectrum template
+name, intensity scale, type of cross dispersion profile, profile
+width in the center line, change of width per line, profile position
+in the center line, and change of position per line.  More specifically:
+
+.ls <template name>
+The one dimensional spectrum template is any one dimensional IRAF image.
+If the spectrum template length is less than the two dimensional spectrum,
+the profile extends only over that number of lines and, if it is longer,
+then only the first part of the spectrum is used.
+.le
+.ls scale
+The template spectrum is scaled by this parameter to define the
+total flux for the two dimensional profile.
+.le
+.ls <profile type>
+The spatial profiles are identified by two keywords, "gaussian"
+or "slit".  The profiles are defined by the following formulae,
+
+.nf
+    gaussian:   I(x) = exp (-ln(2) * (2*(x-xc)/fwhm)**2)
+        slit:   I(x) = exp (-ln(2) * (2*(x-xc)/fwhm)**10)
+.fi
+
+where x is the column coordinate, xc is the profile center, and
+fwhm is the full width at half maximum.  The "gaussian" profile
+is the usual gaussian specified in terms of a FWHM.  The "slit"
+profile is one which is relatively flat and then rapidly drops
+to zero.  The profile is normalized to unit integral so that
+the total flux across the profile is given by the scaled
+1D spectrum flux.
+.le
+.ls fwhm, dfwhm
+The full width at half maximum and derivative with line number.  The fwhm is
+defined for the middle of the image.  The FWHM as a function
+of line, l, is,
+
+	fwhm + (l - nlines/2) * dfwhm
+.le
+.ls center, dcenter
+The profile center and derivative with line number.  The center is
+defined for the middle of the image.  The center as a function
+of line, l, is,
+
+	center + (l - nlines/2) * dcenter
+.le
+
+The provision for having the spectra tilted relative to the columns is
+useful for understanding undersampling effects.  However, note that the
+spectral lines are not perpendicular to the dispersion but are always
+aligned with the image lines.
+.ih
+EXAMPLES
+1. Create an artificial multifiber spectrum:
+
+.nf
+	cl> type multifiber.dat
+	arc 4 gauss 3 0 20 .01
+	spec1 .5 gauss 3 0 30 .01
+	spec2 .4 gauss 3 0 40 .01
+	spec3 .9 gauss 3 0 50 .01
+	spec4 .2 gauss 3 0 60 .01
+	spec5 .6 gauss 3 0 70 .01
+	spec6 1 gauss 3 0 80 .01
+	spec7 1 gauss 3 0 90 .01
+	cl> mk1dspec arc cont=0 peak=500 nl=30
+	cl> mk1dspec spec1 nlines=99 seed=1
+	cl> mk1dspec spec2 nlines=80 seed=2
+	cl> mk1dspec spec3 nlines=45 seed=3
+	cl> mk1dspec spec4 nlines=95 seed=4
+	cl> mk1dspec spec5 nlines=66 seed=5
+	cl> mk1dspec spec6 nlines=90 seed=6
+	cl> mk1dspec spec7 nlines=85 seed=7
+	cl> mk2dspec multifiber model=multifiber.dat
+.fi
+
+In this example artificial one dimensional spectra are generated with
+\fBmk1dspec\fR.
+
+2. Create an artificial multislit spectrum:
+
+.nf
+	cl> type multislit.dat
+	arc 10 slit 18 0 120 .01
+	sky 2.5 slit 18 0 140 .01
+	sky 2.5 slit 18 0 160 .01
+	sky 2.5 slit 18 0 180 .01
+	sky 2.5 slit 18 0 200 .01
+	sky 2.5 slit 18 0 220 .01
+
+	spec1 .05 gauss 3 0 140 .01
+	spec2 .2 gauss 4 0 161 .01
+	spec3 .1 gauss 3 0 179 .01
+	spec4 .1 gauss 3 0 200 .01
+	spec5 .15 gauss 4 0 220 .01
+	cl> mk1dspec sky peak=1 nl=100
+	cl> mk2dspec multislit model=multislit.dat nc=400
+.fi
+
+Note how two spectra are overlaid to provide a sky spectrum with a
+narrower object spectrum.
+
+3. Create an artificial long slit spectrum:
+
+.nf
+	cl> type longslit.dat
+	sky 22 slit 160 0 220 .01 
+	spec5 .05 gauss 3 0 140 .01
+	spec1 .05 gauss 3 0 190 .01
+	spec4 .5 gauss 3 0 220 .01
+	spec2 2 gauss 40 0 220 .01
+	spec5 .1 gauss 3 0 240 .01
+	spec1 .02 gauss 3 0 290 .01
+	cl> mk2dspec longslit model=longslit.dat nc=400
+.fi
+
+Note how objects are overlaid on a long slit sky spectrum.  The width
+of the spec2 spectrum is wider simulating a galaxy spectrum.
+
+4. To include noise use the task \fBmknoise\fR:
+
+.nf
+	cl> mk2dspec longslit model=longslit.dat nc=400
+	cl> mknoise longslit rdnoise=10 gain=2 poisson+ ncos=100
+.fi
+
+5. Use a real long slit spectrum and add an object with an artificial spectrum:
+
+.nf
+	cl> mk1dspec artspec1d nlines=50
+	cl> mk2dspec ls005 out=ls005new model=STDIN
+	artspec1d 1 gauss 5 0 125 0
+	[EOF]
+.fi
+.ih
+SEE ALSO
+mk1dspec, mknoise, mkheader
+.endhelp
diff --git a/noao/artdata/doc/mkechelle.hlp b/noao/artdata/doc/mkechelle.hlp
new file mode 100644
index 00000000..34e022c5
--- /dev/null
+++ b/noao/artdata/doc/mkechelle.hlp
@@ -0,0 +1,585 @@
+.help mkechelle Mar93 noao.artdata
+.ih
+NAME
+mkechelle -- Make artificial echelle spectra
+.ih
+USAGE
+mkechelle images [clobber]
+.ih
+PARAMETERS
+.ls images
+List of echelle spectra to create or modify.
+.le
+.ls clobber (query)
+If an existing image is specified the clobber query parameter is used.
+Normally the parameter is not specified on the command line so that
+a query will be made for each image which exists.  Putting a value
+on the command line permanently overrides the query.  This should be
+done if the task is run in the background.
+.le
+.ls ncols = 512, nlines = 512, norders = 23
+For two dimensional spectra these parameters define the number of columns
+and lines in the final image and the maximum number of orders (there may be
+orders falling outside the image).  The dispersion is along the columns
+which is the second or line axis (dispersion axis is 2) so the number of
+columns is the number of pixels across the dispersion and the number of
+lines is the number of pixels along the dispersion per order.
+
+The extracted format turns the number of lines into the number columns
+and the number of orders is the number of lines; i.e the image
+has the specified number of extracted orders, one per image line,
+with the number of pixels along the dispersion specified by the
+\fInlines\fR parameter.  This is equivalent to what the \fBapextract\fR
+package would  produces for an extracted echelle format with an original
+dispersion axis of 2.  There is no check in this case for orders
+which might fall outside the two dimensional format; i.e. exactly
+the number of orders are created.
+.le
+.ls title = "Artificial Echelle Spectrum"
+Image title to be given to the spectra.  Maximum of 79 characters.
+.le
+.ls header = "artdata$stdheader.dat"
+Image or header keyword data file.  If an image is given then the image
+header is copied.  If a file is given then the FITS format cards are
+copied.  The data file consists of lines in FITS format with leading
+whitespace ignored.  A FITS card must begin with an uppercase/numeric
+keyword.  Lines not beginning with a FITS keyword such as comments or lower
+case are ignored.  The user keyword output of \fBimheader\fR is an
+acceptable data file.  See \fBmkheader\fR for further information.
+.le
+.ls list = no
+List the grating/instrument parameters?
+.le
+.ls make = yes
+Make the artificial spectra?  This is set to no if only the grating
+parameter listing is desired.
+.le
+.ls comments = yes
+Include comments recording task parameters in the image header?
+.le
+
+.ce
+FORMAT PARAMETERS
+.ls xc = INDEF, yc = INDEF
+The column and line position of the blaze peak in the reference order (see
+\fIorder\fR parameter.  If INDEF then the middle of the dimension is used.
+This allows setting the image center relative to the center of the echelle
+pattern.  As with the number of lines and columns the interpretation of
+these numbers relative to the image created depends on whether the format
+is extracted or not.
+.le
+.ls pixsize = 0.027
+Pixel size in millimeters.  This is used to convert the focal length
+and dispersion to pixels.  If INDEF then these parameters are
+assumed to be in pixels.
+.le
+.ls profile = "gaussian" (extracted|gaussian|slit)
+The order profile across the dispersion.  If the value is "extracted"
+then an extracted echelle format spectrum is produced.  Otherwise a
+two dimensional format with a gaussian or slit profile is produced.
+See \fBmk2dspec\fR for a discussion of the profile functions.
+.le
+.ls width = 5.
+If two dimensional echelle images are selected this parameter specifies
+the order profile full width at half maximum in pixels.  See \fBmk2dspec\fR
+for a fuller discussion.
+.le
+.ls scattered = 0.
+Scattered light peak flux per pixel.  A simple scattered light component
+may be included in the two dimensional format.  The scattered light has the
+blaze function shape of the central order along the dispersion and the
+crossdisperser blaze function shape across the dispersion with the peak
+value given by this parameter.  A value of zero indicates no scattered
+light component.
+.le
+
+.ce
+GRATING PARAMETERS
+
+Any of the following parameters may be specified as INDEF.  The missing
+values are resolved using the grating equations described in the
+DESCRIPTION section.  If it is not possible to resolve all the grating
+parameters but the order, wavelength, and dispersion are specified
+then a linear dispersion function is used.  Also in this case the
+extracted format will include dispersion information.
+.ls f = 590., cf = 590.
+Echelle and crossdisperser focal lengths in millimeters (if \fIpixsize\fR
+is given) or pixels.  Technically it is defined by the equation x = f * tan
+(theta) where x is distance from the optical axis on the detector and theta
+is the diffraction angle; i.e. it converts angular measures to millimeters
+or pixels on the detector.  If the focal length is specified as INDEF  it
+may be computed from the dispersion, which is required in this case, and
+the other parameters.
+.le
+.ls gmm = 31.6, cgmm = 226.
+Echelle and crossdisperser grating grooves per millimeter.  If specified as
+INDEF it may be computed from the order, which is required in this case,
+and the other parameters.
+.le
+.ls blaze = 63., cblaze = 4.53
+Echelle and crossdisperser blaze angles in degrees.  It is always specified or printed as a positive
+angle relative to the grating normal.  If specified as INDEF it is
+computed from the other parameters.
+.le
+.ls theta = 69., ctheta = -11.97
+Echelle and crossdisperser angles of incidence in degrees.  The angle of
+incidence must be in the plane perpendicular to face of the grating.  The
+angle of incidence may be specified relative to the grating normal or the
+blaze angle though it is always printed relative to the grating normal.  To
+specify it relative to the blaze angle add 360 degrees; for example to have
+an angle of 15 degrees less than the blaze angle specify 360 - 15 = 345.
+If the angle of incidence is specified as INDEF it is computed from the
+other parameters.
+.le
+.ls order = 112
+The central or reference echelle order for which the wavelength and
+dispersion are specified.  If specified as INDEF it will be computed from
+the grooves per mm, which is required in this case, and the other
+parameters.  In combination with the number of orders this defines the
+first and last orders.  The highest order is the central order plus
+the integer part of one half the number of orders.  However, the
+lowest order is constrained to be at least 1.  The
+reference order is also used in the definitions of \fIxc\fR and \fIyc\fR.
+.le
+.ls corder = 1
+The crossdisperser order for which the crossdisperser blaze wavelength and
+dispersion are specified.  If specified as INDEF it will be computed from
+the grooves per mm, which is required in this case, and the other
+parameters.
+
+If the order is zero then the other grating parameters are ignored and a
+prism-like dispersion is used with the property that the order spacing is
+constant.  Specifically the dispersion varies as the inverse of the
+wavelength with the \fIcwavelength\fR and \fIcdispersion\fR defining the
+function.
+.le
+.ls wavelength = 5007.49 cwavelength = 6700.
+Echelle and crossdisperser blaze wavelengths in Angstroms at the reference
+orders.  If specified as INDEF it will be computed from the other parameters.
+.le
+.ls dispersion = 2.61 cdispersion = 70.
+Echelle and crossdisperser blaze dispersions in Angstroms per millimeter
+(if \fIpixsize\fR is specified) or pixels.
+If specified as INDEF it will be computed from the focal length, which is
+required in this case, and the other parameters.
+.le
+
+.ce
+SPECTRA PARAMETERS
+.ls rv = 0.
+Radial velocity (km/s) or redshift, as selected by the parameter \fIz\fR,
+applied to line positions and continuum.  Velocities are converted to
+redshift using the relativistic relation 1+z = sqrt ((1+rv/c)/(1-rv/c)).
+Note the shift is not a shift in the dispersion parameters but in the
+underlying artificial spectrum.
+.le
+.ls z = no
+Is the velocity parameter a radial velocity or a redshift?
+.le
+.ls continuum = 1000.
+Continuum at the echelle blaze peak in the reference order.
+.le
+.ls temperature = 5700.
+Blackbody continuum temperature in Kelvin.  A value of 0 is used if
+no blackbody continuum is desired.  The intensity level is set by
+scaling to the continuum level at blaze peak reference point.
+.le
+
+.ls lines = ""
+List of spectral line files.  Spectral line files contain lines of rest
+wavelength, peak, and widths (see the DESCRIPTION section).
+The latter two parameters may be missing in which case they default to
+the task \fIpeak\fR and \fIsigma\fR parameters.  If no file or a new
+(nonexistent) file is specified then a number of random lines given by the
+parameter \fInlines\fR is generated.  If a new file name is specified then
+the lines generated are recorded in the file.  If the list of spectral
+line files is shorter than the list of input spectra, the last
+spectral line list file is reused.
+.le
+.ls nlines = 0
+If no spectral line file or a new file is specified then the task will
+generate this number of random spectral lines.  The rest wavelengths are
+uniformly random within the limits of the spectrum, the peaks are
+uniformly random between zero and the value of the \fIpeak\fR parameter
+and the width is fixed at the value of the \fIsigma\fR parameter.
+If a redshift is applied the rest wavelengths are shifted and repeated
+periodically.
+.le
+.ls peak = -0.5
+The maximum spectral line peak value when generating random lines or
+when the peak is missing from the spectral line file.
+This value is relative to the continuum unless the continuum is zero.
+Negative values are absorption lines and positive values are emission lines.
+.le
+.ls sigma = 1.
+The default line width as a gaussian sigma in Angstroms when generating
+random lines or when the width is missing from the spectral line file.
+.le
+.ls seed = 1
+Random number seed.
+.le
+
+PACKAGE PARAMETERS
+.ls nxsub = 10
+Number of pixel subsamples used in computing the gaussian spectral line
+profiles.
+.le
+.ls dynrange = 100000.
+The gaussian line profiles extend to infinity so a dynamic range, the ratio
+of the peak intensity to the cutoff intensity, is imposed to cutoff the
+profiles.
+.le
+.ih
+DESCRIPTION
+This task creates or adds to artificial extracted (one dimensional
+"echelle" format) or two dimensional echelle spectra.  The input spectrum
+(before modification by the spectrograph model) may be a combination of
+doppler shifted blackbody or constant continuum and emission and absorption
+gaussian profile spectral lines.  The lines may have randomly selected
+parameters or be taken from an input file.  Note that the parameters and
+method is similar to the task \fBmk1dspec\fR except that the input line list
+cannot specify a profile type and only Gaussian profiles are currently
+allowed.  The input spectrum is then
+separated out into echelle orders and either recorded as extracted one
+dimensional orders or convolved with a spatial profile and crossdispersed
+into a two dimensional image.  The properties of the echelle grating,
+crossdisperser, and instrumental configuration are modeled described
+later.
+
+If an existing image is specified the \fIclobber\fR parameter is used
+to determine whether to add the generated artificial echelle spectrum
+to the image.  Generally the clobber parameter is not specified on the
+command line to cause a query with the image name to be made for
+each image which already exists.  However, it is possible to put
+the clobber parameter on the command line to eliminate the query.
+This is appropriate for running the task in the background.
+
+There is \fIno\fR checking for consistency with an existing image;
+i.e. that it is an echelle image, whether it is an extracted format
+or a two dimensional spectrum, and what it's wavelength and order
+coverage is.  The only thing that happens is that the \fIncols\fR,
+\fInlines\fR, and \fInorders\fR parameters are replaced by the appropriate
+dimensions of the image with the choice between \fInlines\fR and
+\fInorders\fR made by the \fIprofile\fR parameter (as discussed below)
+and not by the format of the image.
+
+The created spectra are two dimensional, real datatype, images.  A title
+may be given and a set of header keywords be added by specifying an image
+or data file with the \fIheader\fR parameter (see also \fBmkheader\fR).  If
+a data file is specified lines beginning with FITS keywords are entered in
+the image header.  Leading whitespace is ignored and any lines beginning
+with words having lowercase and nonvalid FITS keyword characters are
+ignored.  In addition to this optional header, various parameters which
+occur during reduction of real echelle spectra, such a wavelength
+coordinates for extracted and dispersion corrected spectra, are added.
+Finally, comments may be added to the image header recording the task
+parameters and any information from the line file which are not line
+definitions.
+
+The creation of an artificial echelle spectra has three stages.  First a
+true spectrum is generated; i.e. the spectrum which arrives at the
+spectrograph.  The spectrum is then separated into orders and the
+dispersion and  blaze functions of the echelle and crossdisperser gratings
+(or crossdisperser prism) are applied.  Finally, if a two dimensional
+format is desired it is convolved by an spatial profile (either a gaussian
+or a broader slit-like profile) and the orders are placed as required by
+the crossdispersion relation.
+
+Generation of the model spectrum has three parts; defining a continuum,
+adding emission and absorption lines, and applying a doppler shift.  The
+continuum has two parameters; an intensity scale set by the \fIcontinuum\fR
+parameter and a shape set by the \fItemperature\fR parameter.  The
+intensity scale is set by defining the total, final, extracted intensity in
+a pixel at the blaze peak (rest) wavelength in the reference order; i.e. at
+the wavelength set by the \fIwavelength\fR parameter.  Note this means that
+the efficiency of the gratings at that wavelength is included.  The shape
+of the continuum may be either a blackbody if a positive temperature is
+specified or constant.
+
+Spectral lines are modeled by gaussian profiles of specified wavelength,
+peak, and sigma.  The lines are defined in a spectral line file or
+generated randomly.  A spectral line file consists of text lines giving
+rest wavelength, peak, and sigma.  The sigma or the sigma and peak may be
+absent in which case the parameters \fIsigma\fR and \fIpeak\fR will be
+used.  If peak values are missing random values between zero and the
+\fIpeak\fR value are generated.  Thus, a simple list of wavelengths or a
+list of wavelengths and peaks may be used.
+
+If no spectral line file is specified or a new (nonexistent) file is named
+then the number of random lines given by the parameter \fInlines\fR is
+generated.  The rest wavelengths are uniformly random within the wavelength
+range of the spectrum and extend periodically outside this range in the
+case of an applied velocity shift, the peaks are uniformly random between
+zero and the \fIpeak\fR parameter, and the widths are given by the
+\fIsigma\fR parameter.  If a new file is named then the parameters of the
+generated lines will be output.
+
+The peak values are taken relative to a positive continuum.  In other words
+the generated line profile is multiplied by the continuum (with a minimum
+of zero for fully saturated absorption lines).  If the continuum is less
+than or equal to zero, as in the case of an artificial arc spectrum or pure
+emission line spectrum, then the peak values are interpreted as absolute
+intensities.  Positive peak values produce emission lines and negative
+values produce absorption lines.  Odd results will occur if the continuum
+has both positive and zero or negative values.
+
+The width values are gaussian sigmas given in Angstroms.
+
+The underlying rest spectrum may be shifted.  This is used primarily for
+testing radial velocity measuring algorithms and is not intended as a
+complete model of redshift effects.  The observed wavelength coverage as
+defined by the grating parameters and number of orders is not changed by
+redshifting.  Input line wavelengths are specified at rest and then shifted
+into or out of the final spectrum.  To be realistic the line list should
+include wavelengths over a great enough range to cover all desired
+redshifts.  The peaks and sigma are also appropriately modified by a
+redshift.  As an example, if the redshift is 1 the lines will appear
+broader by a factor of 2 and the peaks will be down by a factor of 2 in
+order to maintain the same flux.
+
+The random line generation is complicated because one wants to have the
+same set of lines (for a given seed) observed at different redshifts.  What
+is done is that the specified number of random lines is generated within
+the observed wavelength interval taken at rest.  This set is then repeated
+periodical over all wavelengths.  A redshift will then shift these rest
+lines in to or out of the observed spectrum.  If the lines are output to a
+line file, they are given at rest.  \fBNote that this periodicity may be
+important in interpreting cross-correlation redshift tests for large shifts
+between template and object spectra.\fR
+
+The definitions of the continuum are also affected by a redshift.  The
+reference point for the continuum level and blackbody shape is the starting
+wavelength taken at rest.  Shifts will then modify the continuum level at
+the reference pixel appropriately.  In particular a large redshift will
+shift the blackbody in such a way that the flux is still given by the
+\fIcontinuum\fR parameter at the reference wavelength at rest.
+
+Once the input spectrum is defined it is modified by the effects of an
+echelle grating and crossdispersion.  This includes the dispersion relation
+between pixel and wavelength, the blaze response function of the gratings,
+and separation into orders.
+
+The primary reference for the model of the echelle grating (a
+crossdisperser grating also obeys this model) used in this task is "Echelle
+efficiencies: theory and experiment" by Schroeder and Hilliard in Applied
+Optics, Vol. 19, No. 16, 1980, p. 2833.  (The nomenclature below is similar
+to that paper except we use theta for alpha, their theta is theta - blaze,
+the reciprocal of the groove spacing which is the grooves per millimeter,
+and the dispersion per linear distance at the detector rather than per
+radian).  This task only treats the case where the incident beam is in the
+plane perpendicular to the grating face (gamma=0).  In this case the basic
+equation is
+
+.nf
+(1)	m * lambda = (sin(theta) + sin(beta)) / g
+.fi
+
+where m is the order, lambda the wavelength, g the grooves per wavelength
+unit, theta the angle of incidence to the grating normal, and beta the
+angle of diffraction to the normal.  The diffraction angle relative to that
+of the blaze maximum, psi, is given by
+
+.nf
+(2)	beta = psi + 2 * blaze - theta
+.fi
+
+where blaze is the blaze angle.  The diffraction angle psi is related to
+position on the detector, again measured from the blaze peak, by
+
+.nf
+(3)	x = f / pixsize * tan(psi)
+.fi
+
+where f is the effective focal length (as defined by this equation) and
+pixsize is the pixel size in millimeters that converts the detector
+positions to pixels.  If a pixel size is not specified then f will be
+taken as being in pixels.
+
+The second basic equation is the diffraction pattern or blaze response
+given by
+
+.nf
+(5)	I = I0 * (sin(delta) / delta) ** 2
+(6)	delta = 2 * pi / lambda * (cos(theta) / g) / cos(epsilon) *
+		sin(psi/2) * cos(epsilon-psi/2)
+(7)	epsilon = theta - blaze
+.fi
+
+where epsilon is the angle between the blaze angle and the angle of
+incidence (the theta of  Shroeder and Hilliard).  When theta = blaze, (6)
+simplifies to
+
+.nf
+(6a)	delta = pi / lambda * (cos (blaze) / g) * sin (psi)
+.fi
+
+As discussed by Schroeder and Hilliard, the relative intensity at the blaze
+peak, I0, must be reduced by the fraction of light at the same wavelength
+as the blaze peak which is diffracted into other orders.  Furthermore at
+some diffraction angles the light is reflected off the second face of the
+grating giving a different effective diffraction angle to be used in (6).
+This computation is done by the task giving a variation in relative blaze
+response with order and reproducing the calculations of Schroeder and
+Hilliard.  The absolute normalization, including the crossdisperser blaze
+function if any, is such that the response at the blaze peak of the
+reference order is unity.  This insures that specified continuum level at
+the reference wavelength is produced.
+
+At the blaze maximum psi = x = 0 and the wavelength and dispersion per
+millimeter on the detector are given by (1) and the derivative of (1) with
+respect to x:
+
+.nf
+(8)	wavelength = 1E7*(sin(theta)+sin(2*blaze-theta))/(gmm*order)
+(9)	dispersion = 1E7*cos(2*blaze-theta)/(gmm*order*f/pixsize)
+.fi
+
+The variable names are the same as the parameters in this task.   In
+particular, gmm is the echelle grooves per millimeter with the factors of
+1E7 (10 to the seventh power) to convert to Angstroms, the factor of f /
+pixsize to convert the dispersion to per pixel, and order is the reference
+order for the wavelength and dispersion.
+
+The \fBmkechelle\fR task provides different ways to define the parameters.
+If there is insufficient information to determine all the grating
+parameters but the wavelength, dispersion, order are specified then
+a simplified grating equation is used which is linear with pixel
+position.  The approximation is that tan(psi) = sin(psi) = psi so
+that
+
+.nf
+(9)	lambda = (order * wavelength + dispersion * x) / m
+               = (a + b * x) / m
+(10)	delta  = pi * order * dispersion / lambda * x
+               =  c / lambda * x
+.fi
+
+Also in this case the extracted format (described later) includes
+wavelength information in the header so that the spectra appear as fully
+dispersion corrected.
+
+If there are at least five of the seven grating parameters specified
+then equations (8) and (9) are used to determine
+unspecified parameters or to override parameters if the equations are
+overspecified.  For example, suppose the grooves per millimeter is known
+but not the blaze angle or focal length.  Then the wavelength and
+dispersion at the reference order are used to compute these quantities.
+
+The full set of grating parameters derived and used to create the spectra
+are documented in the image header if the \fIcomments\fR parameter is
+specified.  Also the \fIlist\fR parameter may be set to print the grating
+parameters and the \fImake\fR parameter may be set to no to check the
+grating parameters without making the spectra.
+
+The crossdisperser grating parameters are treated exactly as above except,
+since only one order is used, the relative blaze efficiency is not
+computed.
+
+There is a variant on the crossdispersion to allow a prism-like separation
+of the echelle orders.  If the crossdispersion grating order, \fIcorder\fR
+is set to zero then the other grating parameters are ignored and a
+prism-like dispersion is used with the property that the order spacing is
+constant.  Specifically the dispersion varies as the inverse of the
+wavelength with the \fIcwavelength\fR and \fIcdispersion\fR defining the
+function.  There is no crossdisperser blaze function in this case either;
+i.e. the relative intensities between orders is solely due to the echelle
+grating blaze response.
+
+There is an interesting effect which follows from the above equations but
+which is not obvious at first glance.  When the full grating equation is
+used the dispersion varies with wavelength.  This means the size of a pixel
+in wavelength varies and so the flux in a pixel changes.  The effect is
+such that the order intensity maximum shifts to the blue from the blaze peak
+because the pixel width in Angstroms increases to the blue faster, for a
+while, than the blaze response decreases.
+
+Once the spectrum has been separated into orders, modified by the
+grating blaze functions, and sampled into pixels in the dispersion
+direction it may be output as an extracted "echelle" format spectrum.
+This occurs when the spatial profile is specified as "extracted".
+The keywords added by the \fBapextract\fR package are included in
+the image header.  If the dispersion model is linear
+the keywords are the same as those produced by the dispersion
+correction task \fBecdispcor\fR.
+
+If the spatial profile is specified as "gaussian" or "slit" then the
+orders are convolved by the profile function and the crossdispersion
+relation is used to determine where the order falls at each wavelength.
+The spatial profiles are defined by the formulas:
+
+.nf
+    gaussian:   I(x) = exp (-ln(2) * (2*(x-xc(w))/width)**2)
+        slit:   I(x) = exp (-ln(2) * (2*(x-xc(w))/width)**10)
+.fi
+
+where x is the spatial coordinate, xc(w) is the order center at
+wavelength w, and width is the full width at half maximum specified by
+the parameter of that name.  The "gaussian" profile
+is the usual gaussian specified in terms of a FWHM.  The "slit"
+profile is one which is relatively flat and then rapidly drops
+to zero.  The profile is normalized to unit integral so that
+the total flux across the profile is given by the scaled
+1D spectrum flux.  The profile is fully sampled and then binned to
+the pixel size to correctly include sampling effects as a function
+of where in a pixel the order center falls.
+
+Note that in this model the orders are always tilted with respect
+to the columns and constant wavelength is exactly aligned with the
+image lines.
+.ih
+EXAMPLES
+1. Create an absorption spectrum with blackbody continuum and scattered
+light using the default grating parameters then add noise.
+
+.nf
+	cl> mkechelle ex1 nrand=100 scat=100.
+	cl> mknoise ex1 gain=2 rdnoise=5 poisson+
+.fi
+
+2. Create an arc spectrum using the line list noao$lib/onedstds/thorium.dat.
+
+.nf
+	cl> mkechelle ex2 cont=10 temp=0 \
+	lines=noao$lib/onedstds/thorium.dat peak=1000 sigma=.05
+.fi
+
+Note that the line intensities are random and not realistic.  The peak
+intensities range from 0 to 1000 times the continuum or 10000.
+
+3. Create an extracted version of example1.
+
+.nf
+	cl> mkechelle ex1.ec prof=extracted nrand=100 scat=100.
+	cl> mknoise ex1.ec gain=2 rdnoise=5 poisson+
+.fi
+
+Note that the noise is different and greater than would be the case with
+extracting the orders of example 1 because the noise is not summed
+over the order profile but is added after the fact.
+
+4. Create an extracted and dispersion corrected version of example1.
+
+.nf
+	cl> mkechelle ex1.ec prof=extracted nrand=100 scat=100. \
+	gmm=INDEF blaze=INDEF theta=INDEF
+	Echelle grating: Using linear dispersion
+	Warning: Insufficient information to resolve grating parameters
+	cl> mknoise ex1.ec gain=2 rdnoise=5 poisson+
+.fi
+
+The warning is expected.  By not specifying all the parameters needed to
+fully model an echelle grating the default action is to use a linear
+dispersion in each order and to set the image header dispersion
+information.  When a complete grating model is specified, as in example 3,
+the extracted spectrum is not given dispersion information so that the
+nonlinear behavior of the dispersion can be applied by \fBecidentify\fR and
+\fBdispcor\fR.  As with example 3, the noise is different since it is added
+after extraction and dispersion correction.
+.ih
+REVISIONS
+.ls MKECHELLE V2.10.3
+The task was updated to produce the current coordinate system format.
+.le
+.ih
+SEE ALSO mknoise, mk1dspec, mk2dspec, mkheader, astutil.gratings
+.endhelp
diff --git a/noao/artdata/doc/mkexamples.hlp b/noao/artdata/doc/mkexamples.hlp
new file mode 100644
index 00000000..b3b012a0
--- /dev/null
+++ b/noao/artdata/doc/mkexamples.hlp
@@ -0,0 +1,167 @@
+.help mkexamples Mar93 noao.artdata
+.ih
+NAME
+mkexamples - Make artificial data example images
+.ih
+USAGE
+mkexamples name image
+.ih
+PARAMETERS
+.ls name
+Example name (abbreviations are not allowed):
+
+.nf
+     galcluster - Galaxy cluster
+       globular - Globular cluster
+       galfield - Galaxy field
+      starfield - Starfield
+
+         henear - Helium-neon-argon spectrum (uncalibrated)
+       heneardc - Helium-neon-argon spectrum (calibrated)
+
+          ecarc - Echelle thorium-argon spectrum (uncalibrated)
+        ecarcdc - Echelle thorium-argon spectrum (calibrated)
+
+       spectrum - Absorption spectrum (calibrated)
+        echelle - Echelle absorption spectrum (calibrated)
+
+        ecarc2d - Echelle thorium-argon slit spectrum
+        ecobj2d - Echelle object slit spectrum
+
+          lsarc - Long slit helium-neon-argon spectrum
+	  lsobj - Long slit object spectrum
+
+     multifiber - Multifiber spectrum
+.fi
+.le
+.ls image
+Output image name.
+.le
+.ls oseed = 1
+Random number seed affecting object generation.  Different object seeds
+will produces different examples of objects or spectral lines or number
+of apertures/orders.  This
+usually modifies the seed parameters in \fBstarlist\fR, \fBgallist\fR,
+and \fBmk1dspec\fR.
+.le
+.ls nseed = 1
+Random number noise seed.  Different noise seeds will produce examples
+with different noise, generally of the same level but simply having
+a different pattern.  This is usually the seed parameter in
+\fBmkobjects\fR or \fBmknoise\fR.
+.le
+.ls comments = no
+Add comments to the image header describing various artificial data
+parameters?
+.le
+.ls verbose = yes
+Print message indicating image being created?
+.le
+.ls errors = yes
+Print messages if the image already exists, bad example name, or other
+errors?
+.le
+.ls list = no
+List script used to generate the example rather than create an image?
+.le
+.ih
+DESCRIPTION
+The task is intended to generate a few artificial images of various types to
+be used as examples of the artificial data package and in various
+demonstrations and test procedures for other packages.  The examples are not 
+exhaustive.  The only adjustable parameters are variations of the
+random number seeds.  Varying the noise seed allows several observations
+of the same example while varying the object seed allows several observations
+of different "fields", spectral lines, or number of apertures/orders.
+
+If the example name is not given on the command line a menu of example
+names is first printed and then a prompt for the name is given.
+The name may be a submenu or an example.  The
+names may not be abbreviated.  If desired the simple command
+script used to generate the example may be paged.  Otherwise the
+specified image will be generated.  Keep in mind that some of the
+examples (particularly those generating galaxy images) may take a
+significant amount of time.  On a SPARCstation the examples all run in
+under five minutes.  A check is made to see if the image already
+exists.  If the image exists then the task exits.  If the \fIerrors\fR
+parameter is specified an error message is printed.
+
+A reason for the error output to be turned off is in test scripts and
+demonstrations where the image will be created the first time and reused
+in further tests or demonstrations.  In such cases the verbose option is
+generally set so that the user is aware that an image is being created
+and some delay is to be expected.
+
+This task is a procedure script which selects and lists or executes
+any file in the mkexamples$ logical directory with the example name and the
+extension ".cl".  Thus, to add additional examples create a simple
+command script (not a procedure script) and place it in the mkexamples
+directory along with an entry in the menu file mkexamples$mkexamples.men.
+.ih
+EXAMPLES
+1. Create a globular cluster example.
+
+.nf
+    ar> mkexample
+	    		MKEXAMPLE Menu
+
+     galcluster - Galaxy cluster
+       globular - Globular cluster
+       galfield - Galaxy field
+      starfield - Starfield
+
+       onedspec - Menu of one dimensional spectra
+       twodspec - Menu of two dimensional spectra
+     threedspec - Menu of three dimensional spectra
+    Example name: globular
+    Image name: globular
+    Creating example globular in image globular ...
+.fi
+
+2.  Try and create the same example again.
+
+.nf
+    ar> mkexample globular globular
+    ERROR: Image globular already exists
+.fi
+
+3.  List the script which creates the globular example.
+
+.nf
+    ar> mkexample globular list+
+    # GLOBULAR - Globular cluster
+
+    file	image, dat
+
+    image = s1
+    dat = mktemp ("art")
+
+    starlist (dat, 5000, "", "", interactive=no, spatial="hubble",
+	xmin=1., xmax=512., ymin=1., ymax=512., xcenter=INDEF,
+	ycenter=INDEF, core_radius=30., base=0., sseed=i,
+	luminosity="bands", minmag=-7., maxmag=0., mzero=-4., power=0.6,
+	alpha=0.74, beta=0.04, delta=0.294, mstar=1.28, lseed=i,
+	nssample=100, sorder=10, nlsample=100, lorder=10,
+	rbinsize=10., mbinsize=0.5, graphics="stdgraph", cursor="")
+
+    mkobjects (image, output="", ncols=512, nlines=512,
+	title="Example artificial globular cluster",
+	header="artdata$stdheader.dat", background=1000., objects=dat,
+	xoffset=0., yoffset=0., star="moffat", radius=1.0, beta=2.5,
+	ar=1., pa=0., distance=1., exptime=1., magzero=7.,
+	gain=3., rdnoise=10., poisson=yes, seed=j)
+
+    delete (dat, verify=no)
+.fi
+.ih
+REVISIONS
+.ls MKEXAMPLES V2.10.3
+The examples have been expanded to include submenus.  The submenus organize
+the various types of spectra.  Additional spectral examples have been
+added.  The oseed parameter selects the number of apertures in the
+onedspec spectra and the number of orders in the echelle examples.
+.le
+.ih
+SEE ALSO
+mkobjects, mknoise, mk1dspec, mk2dspec, mkechelle
+.endhelp
diff --git a/noao/artdata/doc/mkheader.hlp b/noao/artdata/doc/mkheader.hlp
new file mode 100644
index 00000000..43b37895
--- /dev/null
+++ b/noao/artdata/doc/mkheader.hlp
@@ -0,0 +1,84 @@
+.help mkheader Aug90 noao.artdata
+.ih
+NAME
+mkheader - Append/replace image header
+.ih
+USAGE
+mkheader images headers
+.ih
+PARAMETERS
+.ls images
+List of images in which header information is to be added or modified.
+.le
+.ls header = "artdata$stdheader.dat"
+List of images or header keyword data files.  If the list is shorter
+than the input image list then the last entry is repeated.
+If an image is given then the image header
+is copied.  If a file is given then the FITS format cards are copied.
+This only applies to new images.   The data file consists of lines
+in FITS format with leading whitespace ignored.  A FITS card must begin
+with an uppercase/numeric keyword.  Lines not beginning with a FITS
+keyword such as comments or lower case are ignored.  The user keyword
+output of \fBimheader\fR is an acceptable data file.
+.le
+.ls append = yes
+Append to existing keywords?  If no then the existing header is replaced.
+.le
+.ls verbose = no
+Verbose output?
+.le
+.ih
+DESCRIPTION
+The image headers in the list of input images may be replaced or appended
+with information from images or data files specified by the \fIheader\fR
+parameter list.  If the header list is shorter than the list of images
+to be modified the last header file is repeated.  Depending on the
+value of the \fIappend\fR parameter, new parameters will be appended
+or replace the existing image header parameters.
+
+A header keyword data file consists of lines of FITS format cards.
+Leading whitespace is ignored.  Lines not recognized as FITS cards
+are ignored.  A valid FITS card is defined as beginning with a keyword
+of up to 8 uppercase, digit, hyphen, or underscore characters.  If
+less than 8 characters the remaining characters are blanks.  The
+ninth character may be an equal sign but must be immediately followed
+by a blank.  Such value cards should be in FITS format though no
+attempt is made to enforce this.  Any other ninth character is also
+acceptable and the line will be treated as a comment.  Note that this
+way of recognizing FITS parameters excludes the case of comments
+in which the first 8 characters are blank.  The reason for allowing
+leading whitespace and eliminating the blank keyword case is so that
+the long output of \fBimheader\fR may be used directly as input.
+
+Header files are also used by several of the tasks in the artificial
+data package with a standard default file "artdata$stdheader.dat".
+To edit image headers also see \fBhedit\fR.
+.ih
+EXAMPLES
+1. Add some standard keywords from a file to an image.
+
+.nf
+    ar> type myheader
+    # MY header list
+    INSTRUME= 'bspec mark II'		/ B Spectrograph
+    LENS    =                  3	/ Lens number
+    FOCRATIO=                5.2        / Focal ratio
+    ar> mkheader *.imh myheader
+.fi
+
+2. Copy an image header.
+
+    ar> mkheader new dev$pix append-
+
+3. Edit the image header with a text editor and replace the old header
+with the edited header.
+
+.nf
+    ar> imheader myimage l+ > temp
+    ar> edit temp
+    ar> mkheader myimage temp append-
+.fi
+.ih
+SEE ALSO
+hedit, mkobjects, mknoise, mk1dspec, mk2dspec
+.endhelp
diff --git a/noao/artdata/doc/mknoise.hlp b/noao/artdata/doc/mknoise.hlp
new file mode 100644
index 00000000..cab14d42
--- /dev/null
+++ b/noao/artdata/doc/mknoise.hlp
@@ -0,0 +1,245 @@
+.help mknoise Aug90 noao.artdata
+.ih
+NAME
+mknoise - Make/add noise and cosmic rays to 1D/2D images
+.ih
+PARAMETERS
+.ls input
+Images to create or modify.
+.le
+.ls output = ""
+Output images when modifying input images.  If no output images are
+given then existing images in the input list are modified directly.
+If an output image list is given then it must match in number the
+input list.
+.le
+
+WHEN CREATING NEW IMAGES
+.ls title = ""
+Image title to be given to the images.  Maximum of 79 characters.
+.le
+.ls ncols = 512, nlines = 512
+Number of columns and lines.
+.le
+.ls header = "artdata$stdheader.dat"
+Image or header keyword data file.  If an image is given then the image header
+is copied.  If a file is given then the FITS format cards are copied.
+This only applies to new images.   The data file consists of lines
+in FITS format with leading whitespace ignored.  A FITS card must begin
+with an uppercase/numeric keyword.  Lines not beginning with a FITS
+keyword such as comments or lower case are ignored.  The user keyword
+output of \fBimheader\fR is an acceptable data file.  See \fBmkheader\fR
+for further information.
+.le
+
+NOISE PARAMETERS
+.ls background = 0.
+Background to add to images before computing Poisson noise.
+.le
+.ls gain = 1.
+Gain in electrons per data number.  The gain is used for scaling the
+read noise parameter and in computing poisson noise.
+.le
+.ls rdnoise = 0.
+Gaussian read noise in electrons.
+.le
+.ls poisson = no
+Add poisson noise?  Note that any specified background is added to new
+or existing images before computing the Poisson noise.
+.le
+.ls seed = 1
+Random number seed.  If a value of "INDEF" is given then the clock
+time (integer seconds since 1980) is used as the seed yielding
+different random numbers for each execution.
+.le
+
+COSMIC RAYS
+.ls cosrays = ""
+List of cosmic ray files.  Cosmic ray files contain lines of cosmic ray
+coordinates and energy (see DESCRIPTION section).  If no
+file or a new (nonexistent) file is specified then a number of random
+cosmic rays given by the parameter \fIncosrays\fR is generated.  If a
+new file name is specified then the events generated are recorded in the
+file.  If the list of cosmic ray files is shorter than the list of
+input images then the last cosmic ray file is reused.
+.le
+.ls ncosrays = 0
+If no cosmic ray file or a new file is specified then the task will
+generate this number of random cosmic rays.  The positions are
+uniformly random within the limits of the image and the energy is
+uniformly random between zero and a maximum.
+.le
+.ls energy = 30000.
+When generating random events the cosmic rays will have a uniform energy
+distribution (in electrons) between zero and this maximum.
+.le
+.ls radius = 0.5
+The half-intensity radius of gaussian profile cosmic rays in pixels
+along the major axis.
+.le
+.ls ar = 1.
+Minor to major axial ratio for cosmic rays.
+.le
+.ls pa = 0.
+Position angle in degrees measured counterclockwise from the X axis for
+cosmic rays.
+.le
+
+.ls comments = yes
+Include comments recording task parameters in the image header?
+.le
+
+PACKAGE PARAMETERS
+
+These parameters define certain computational shortcuts which greatly
+affect the computational speed.  They should be adjusted with care.
+.ls nxc = 5, nyc = 5
+Number of cosmic ray centers per pixel in X and Y.  Rather than evaluate
+cosmic rays precisely at each subpixel coordinate, a set of templates
+with a grid of subpixel centers is computed and then the nearest template to
+the desired position is chosen.  The larger the number the more memory
+and startup time required.
+.le
+.ls nxsub = 10, nysub = 10
+Number of pixel subsamples in X and Y used in computing the cosmic
+ray profiles.  This is the subsampling in the central
+pixel and the number of subsamples decreases linearly from the center.
+This affects the time required to compute the cosmic ray templates.
+.le
+.ls dynrange = 100000.
+The intensity profile of the gaussian cosmic rays extends to infinity so
+a dynamic range, the ratio of the peak intensity to the cutoff
+intensity, is imposed.  Because the cosmic rays are small this parameter
+is not critical.
+.le
+.ls ranbuf = 0
+Random number buffer size.  When generating readout and poisson noise,
+evaluation of new random values has an affect on the execution time.
+If truly (or computationally truly) random numbers are not needed
+then this number of random values is stored and a simple
+uniform random number is used to select from the stored values.
+To force evaluation of new random values for every pixel set the
+value of this parameter to zero.
+.le
+.ih
+DESCRIPTION
+This task creates or modifies images with readout noise, poisson noise,
+and cosmic ray events.  New images are created with the specified
+dimensions and real datatype.  Existing images may be modified in place
+or new images may be created.
+
+If a new image is created it is has the mean level given by the parameter
+\fIbackground\fR.  With no noise and no cosmic rays this task can be used to
+create images of constant background value.  For existing images the
+background is added before computing any noise.  To add noise to an
+existing image without modifying the mean counts set the background
+to zero.
+
+For new images a set of header keywords may be added by specifying an
+image or data file with the \fIheader\fR parameter (see also \fBmkheader\fR).
+If a data file is specified lines beginning with FITS keywords are
+entered in the image header.  Leading whitespace is ignored and any
+lines beginning with words having lowercase and nonvalid FITS keyword
+characters are ignored.  In addition to this optional header,
+keywords, parameters for the gain and read noise are defined.
+Finally, comments may be added to the image header recording the task
+parameters and any information from the cosmic ray file which are not
+cosmic ray definitions.
+
+Poisson photon noise is generated by setting the \fIpoisson\fR parameter.
+For new images the input data value is the background while for
+existing images the input data value is added to the background value.
+The data value is then multiplied by the gain, a poisson deviate is
+generated, and divided by the gain.  Expressed as a formula:
+
+.nf
+      New images: out = P(background * gain) / gain
+ Existing images: out = P((in+background)*gain) / gain
+.fi
+
+where P(x) is a poisson deviate with mean x, in and out are the input
+and final pixel values, and background and gain are the parameter
+values of the same name.
+
+Readout or gaussian noise is generated by specifying a gaussian sigma with
+the parameter \fIrdnoise\fR.  The sigma is divided by the specified gain
+to convert to image data units.  Gaussian random numbers of mean zero are
+then generated for each pixel and added to the image, or background
+value for new images, after the photon noise is computed.
+
+Generating gaussian and poisson random numbers computationally is
+the main determinant of the execution time in this task.
+Two things are done to speed up the task.
+First, the gaussian approximation is used for data values greater
+than 20 (after applying the background and gain).  The square root
+of the data value is used as the gaussian sigma about the data
+value.  For values less than 20 a true poisson deviate is generated.
+The second speed up is to allow storing a number of normalized gaussian
+values given by the package parameter \fIranbuf\fR as they are generated.  If
+more values than this are desired then a uniform random number is used
+to select one of these stored values.  This applies to both the read noise
+and poisson noise gaussian approximation though not the true poisson
+evaluation.  For most purposes this approximation is good and one would
+need to look very hard to detect the nonrandomness in the noise.
+However, if one wants to take the extra computational time then
+by setting the \fIranbuf\fR parameter to zero each gaussian
+random number will be generated independently.
+
+The cosmic ray model is an elliptical gaussian of specified
+half-intensity radius, axial ratio, and position angle.  Normally the
+radius will be small (smaller than the point spread function) and the
+axial ratio will be 1.  The cosmic rays are subsampled and can have the
+number of centers given by the \fInxc/nyc\fR package parameters.  The method
+of generating the cosmic rays is that described for the task
+\fBmkobjects\fR.  Specifically it is the same as adding gaussian
+profile stars.
+
+The total flux (not the peak) of the cosmic ray is given by the energy
+in electrons so that the value is divided by the gain to produce the
+total flux in the image.  Note that this task can be used to add cosmic
+ray spikes to one dimensional images such as spectra but the strengths
+will appear low because of the part of the event which falls outside
+the single line.
+
+The positions and energies of the cosmic rays can be specified in a
+file or the task can generate random events.  Specific cosmic rays are
+specified by a file containing lines of x and y positions and energy.
+Positions outside the limits of the image are ignored.  If no cosmic
+ray file is given or if a new, nonexistent file is named then the
+number of cosmic rays given by the \fIncosrays\fR parameter is
+generated with uniform spatial distribution within the image and
+uniform energy distribution between zero and that given by the
+\fIenergy\fR parameter.  By giving a new file name the randomly
+generated cosmic rays will be recorded for reuse or to allow
+identifying the events while testing tasks and algorithms.
+.ih
+EXAMPLES
+1. Create a new image with a background of 1000, a read noise
+of 10 electrons, a gain of 2, and 50 random cosmic rays.  Don't keep a
+record of the cosmic rays.
+
+.nf
+	cl> mknoise testim back=1000 rd=10 gain=2 poisson+ ncos=50
+.fi
+
+2. Add cosmic rays to an image and create a new output image.
+
+.nf
+	cl> head cosfile
+	20.3 50.1 1000
+	325.6 99.6 250
+	cl> mknoise dev$pix out=newpix cos=cosfile
+.fi
+.ih
+REVISIONS
+.ls MKNOISE V2.11+
+The random number seed can be set from the clock time by using the value
+"INDEF" to yield different random numbers for each execution.
+.le
+.ls MKNOISE V2.11
+The default value of "ranbuf" was changed to zero.
+.le
+.ih
+SEE ALSO
+mkobjects, mkheader
+.endhelp
diff --git a/noao/artdata/doc/mkobjects.hlp b/noao/artdata/doc/mkobjects.hlp
new file mode 100644
index 00000000..e641635a
--- /dev/null
+++ b/noao/artdata/doc/mkobjects.hlp
@@ -0,0 +1,636 @@
+.help mkobjects Jan92 noao.artdata
+.ih
+NAME
+mkobjects - Make/add artificial stars and galaxies to 2D images
+.ih
+USAGE
+mkobjects input
+.ih
+PARAMETERS
+.ls input
+Images to create or modify.
+.le
+.ls output = ""
+Output images when modifying input images.  If no output images are
+given then existing images in the input list are modified directly.
+If an output image list is given then it must match in number the
+input list.
+.le
+
+WHEN CREATING NEW IMAGES
+.ls title = ""
+Image title to be given to the images.  Maximum of 79 characters.
+.le
+.ls ncols = 512, nlines = 512
+Number of columns and lines.
+.le
+.ls header = "artdata$stdheader.dat"
+Image or header keyword data file.  If an image is given then the image header
+is copied.  If a file is given then the FITS format cards are copied.
+This only applies to new images.   The data file consists of lines
+in FITS format with leading whitespace ignored.  A FITS card must begin
+with an uppercase/numeric keyword.  Lines not beginning with a FITS
+keyword such as comments or lower case are ignored.  The user keyword
+output of \fBimheader\fR is an acceptable data file.  See \fBmkheader\fR
+for further information.
+.le
+.ls background = 1000.
+Default background and poisson noise background.  This is in data numbers
+with the conversion to photons determined by the \fIgain\fR parameter.
+.le
+
+OBJECT PARAMETERS
+.ls objects = ""
+List of object files.  The number of object files must match the number of
+input images.  The object files contain lines of object coordinates,
+magnitudes, and shape parameters (see the DESCRIPTION section).
+.le
+.ls xoffset = 0.,  yoffset = 0.
+X and Y coordinate offset to be added to the object list coordinates.
+.le
+.ls star = "moffat"
+Type of star and point spread function.  The choices are:
+.ls gaussian
+An elliptical Gaussian profile with major axis half-intensity radius
+given by the parameter \fIradius\fR, axial ratio given by the parameter
+\fIar\fR, and position angle given by the parameter \fIpa\fR.
+.le
+.ls moffat
+An elliptical Moffat profile with major axis half-intensity radius
+given by the parameter \fIradius\fR, model parameter \fIbeta\fR,
+axial ratio given by the parameter \fIar\fR, and position angle given
+by the parameter \fIpa\fR.
+.le
+.ls <image>
+If not one of the profiles above, an image of the specified name is
+sought.  If found the center of the template image is assumed to be the
+center of the star/psf and the image template is scaled so that the
+radius of the template along the first axis is given by the \fIradius\fR
+parameter.  The axial ratio and position angle define an
+elliptical sampling of the template.
+.le
+.ls <profile file>
+If not one of the above, a text file is sought giving either an intensity
+per unit area profile or a cumulative flux profile from the center to the
+edge.  The two are differentiated by whether the first profile point is 0
+for a cumulative profile or nonzero for an intensity profile.  An intensity
+profile is recommended.  If found the profile defines an elliptical star/psf
+with the major axis radius to the last profile point given by the parameter
+\fIradius\fR, axial ratio given by the parameter \fIar\fR, and position
+angle given by the parameter \fIpa\fR.
+.le
+.le
+.ls radius = 1.
+Seeing radius/scale in pixels along the major axis.  For the "gaussian"
+and "moffat" profiles this is the half-intensity radius of the major
+axis, for image templates this is the template radius along the x dimension,
+specifically one half the number of columns, and for arbitrary user profiles
+this is the radius to the last profile point.
+.le
+.ls beta = 2.5
+Moffat model parameter.  See the DESCRIPTION for a definition of the
+Moffat profile.
+.le
+.ls ar = 1.
+Minor to major axial ratio for the star/psf.
+.le
+.ls pa = 0.
+Position angle in degrees measured counterclockwise from the X axis
+for the star/psf.
+.le
+.ls distance = 1.
+Relative distance to be applied to the object list coordinates,
+magnitudes, and scale sizes.  This factor is divided into the
+object coordinates, after adding the offset factors, to allow expanding
+or contracting about any origin.  The magnitudes scale as the
+square of the distance and the sizes of the galaxies scale
+linearly.  This parameter allows changing image sizes and fluxes
+at a given seeing and sampling with one value.
+.le
+.ls exptime = 1.
+Relative exposure time.  The object magnitudes and background
+level are scaled by this parameter.  This is comparable to changing the
+magnitude zero point except that it includes changing the background.
+.le
+.ls magzero = 7.
+Magnitude zero point defining the conversion from magnitudes in the
+object list to instrumental/image fluxes.
+.le
+
+NOISE PARAMETERS
+.ls gain = 1.
+Gain in electrons per data number.  The gain is used for scaling the
+read noise parameter, the background, and in computing poisson noise.
+.le
+.ls rdnoise = 0.
+Gaussian read noise in electrons.  For new images this applies to the
+entire image while for existing images this is added only to the objects.
+.le
+.ls poisson = no
+Add poisson photon noise?  For new images this applies to the entire image
+while for existing images this is only applied to the objects.  Note
+that in the latter case the background parameter is added before
+computing the new value and then subtracted again.
+.le
+.ls seed = 1
+Random number seed.  If a value of "INDEF" is given then the clock
+time (integer seconds since 1980) is used as the seed yielding
+different random numbers for each execution.
+.le
+
+.ls comments = yes
+Include comments recording task parameters in the image header?
+.le
+
+PACKAGE PARAMETERS
+
+These parameters define certain computational shortcuts which greatly
+affect the computational speed.  They should be adjusted with care.
+.ls nxc = 5, nyc = 5
+Number of star and psf centers per pixel in X and Y.  Rather than evaluate
+stars and the psf convolution functions precisely at each subpixel
+coordinate, a set of templates with a grid of subpixel centers is
+computed and then the nearest template to the desired position is chosen.
+The larger the number the more memory and startup time required.
+.le
+.ls nxsub = 10, nysub = 10
+Number of pixel subsamples in X and Y used in computing the star and
+psf.  This is the subsampling in the central
+pixel and the number of subsamples decreases linearly from the center.
+The larger the numbers the longer it takes to compute the star and psf
+convolution templates.
+.le
+.ls nxgsub = 5, nygsub = 5
+Number of pixel subsamples in X and Y used in computing galaxy images.
+This is the subsampling in the central pixel and the number of
+subsamples decreases linearly from the center.  Because galaxy images
+are extended and each subsample is convolved by the psf convolution it
+need not be as finely sampled as the stars.  This is a critical
+parameter in the execution time if galaxies are being modeled.
+The larger the numbers the longer the execution time.
+.le
+.ls dynrange = 100000., psfrange = 10.
+The intensity profiles of the analytic functions extend to infinity so
+a dynamic range, the ratio of the peak intensity to the cutoff
+intensity, is imposed to cutoff the profiles.  The \fIdynrange\fR
+parameter applies to the stellar templates and to the galaxy profiles.
+The larger this parameter the further the profile extends.
+When modeling galaxies this has a fairly
+strong affect on the time (larger numbers means larger images and more
+execution time).  Only for very high signal-to-noise
+objects will the cutoff be noticeable.  A correction is made to
+the object magnitudes to reflect light lost by this cutoff.
+
+The psf convolution, used on galaxies, is generally not
+evaluated over as large a dynamic range, given by the parameter
+\fIpsfrange\fR, especially since it has a very strong affect on the
+execution time.  The convolution is normalized to unit weight over the
+specified dynamic range.
+.le
+.ls ranbuf = 0
+Random number buffer size.  When generating readout and poisson noise,
+evaluation of new random values has an affect on the execution time.
+If truly (or computationally truly) random numbers are not needed
+then this number of random values is stored and a simple
+uniform random number is used to select from the stored values.
+To force evaluation of new random values for every pixel set the
+value of this parameter to zero.
+.le
+.ih
+DESCRIPTION
+This task creates or modifies images by adding models of astronomical
+objects, stars and galaxies, as specified in object lists.  New images are
+created with the specified dimensions, background, title, and real datatype.
+Existing images may be modified in place or new images output.  The
+task includes the effects of image scale, pixel sampling, atmospheric
+seeing, and noise.  The object models may be analytic one dimensional
+profiles, user defined one dimensional profiles, and user defined image
+templates.  The profiles and templates are given elliptical shapes by
+specifying a scale radius for the major axis, a minor axis to major
+axis axial ratio, and a position angle.
+
+For new images a set of header keywords may be added by specifying an
+image or data file with the \fIheader\fR parameter (see also \fBmkheader\fR).
+If a data file is specified lines beginning with FITS keywords are
+entered in the image header.  Leading whitespace is ignored and any
+lines beginning with words having lowercase and nonvalid FITS keyword
+characters are ignored.  In addition to this optional header,
+keywords, parameters for the gain, read noise, and exposure time are
+defined.  Finally, comments may be added to the image header recording the task
+parameters and any information from the objects file which are not
+object definitions; in particular, the \fBstarlist\fR and
+\fBgallist\fR parameters are recorded.
+
+A completely accurate simulation of the effects of pixel sampling,
+atmospheric seeing, object appearance, luminosity functions, and noise
+can require a large amount of computer time even on
+supercomputers.  This task is intended to allow generation of large
+numbers of objects and images over large image sizes representative of
+current deep optical astronomical images.  All this is to be done
+on typical workstations.  Thus, there are many approximations and
+subtle algorithms used to make this possible to as high a degree of
+accuracy as practical.  The discussion will try to describe these in
+sufficient detail for the user to judge the accuracy of the artificial
+data generated and understand the trade offs with many of the
+parameters.
+
+New images are created with the specified dimensions, title, and real
+datatype.  The images have a constant background value given by the
+\fIbackground\fR parameter (in data numbers) before adding objects and
+noise.  Noise consists of gaussian and poisson components.  For existing
+images, noise is only added to the objects and the background parameter is
+used in the calculation of the poisson noise: specifically, a poisson
+random value with mean given by the sum of the object and the background is
+generated and then the background is subtracted.  For more on how the noise
+is computed and approximations used see \fBmknoise\fR.
+
+Objects are specified by a position, magnitude, model, scale, axial
+ratio, and position angle.  Since the point spread function (PSF)
+is assumed constant over the image the star model, size, axial ratio,
+and position angle are specified by the task parameters \fIstar\fR,
+\fIradius\fR, \fIar\fR, and \fIpa\fR.  For galaxies, where the
+intrinsic shapes vary from object to object, these parameters are
+specified as part of the object lists.  For both types of objects the
+positions and magnitudes are specified in the object lists.
+
+There is a great deal of flexibility in defining the object models.
+The models are defined either in terms of a one dimensional radial
+intensity or cumulative flux profile
+or an image template.  The flux profiles may be
+analytic functions or a user defined profile given as an equally spaced
+set of values in a text file.  The first point is zero at the center
+for a cumulative profile
+and increases monotonically to the edge.  Note that intensity profiles
+are to be preferred to avoid artifacts in the conversion from cumulative
+flux.  In particular, cumulative flux profiles may give a spike at the
+center.  In either case, the profile should be specified fairly finely,
+many points, to avoid interpolation effects.
+
+The functional form of the analytic profiles the user profiles, and
+image template are given below.
+
+.nf
+      gaussian:  I = exp (-ln (2) * (R/radius)**2)
+        moffat:  I = (1 + (2**(1/beta)-1) * (R/radius)**2) ** -beta
+     sersic<n>:  I = exp (-b * (R/radius)**1/n)
+       expdisk:  I = exp (-1.6783 * R/radius)
+        devauc:  I = exp (-7.67 * (R/radius)**1/4)
+  flux profile:  I = intensity (nprofile * R/radius)
+  flux profile:  F = flux (nprofile * R/radius)
+image template:  I = image (nc/2+nc/2*dX/radius, nl/2+nc/2*dY/radius)
+.fi
+
+where R, dX, and dY are defined below, \fIradius\fR is the scale parameter
+and \fIbeta\fR is the Moffat parameter specified by the user,
+nprofile is the number of profile points in the user profile, and nc and nl
+are the image template column and line dimensions.  The Gaussian, "gaussian",
+and Moffat, "moffat", profiles are used for stars and the point spread
+function, while the Sersic (sersic),  exponential disk (expdisk), and
+De Vaucouleurs (devauc) profiles are common models for spiral and elliptical
+galaxies.  The image templates are intended to model images with
+some complex structure.  The usual case is to have a very well sampled
+and high signal-to-noise image be reduced in scale (a more distant
+example), convolved with seeing (loss of detail), and noise (degraded
+signal-to-noise).  This also allows for more complex point spread
+functions.
+
+The radial profiles are mapped into two dimensional objects by an elliptical
+transformation.  The image templates are also mapped by an elliptical
+transformation to rotate and stretch them.  If the output image
+coordinates are given by (x, y), and the specified object center
+coordinates are given by (xc, yc) then the transformation is defined
+as shown below.
+
+.nf
+	dx = x - xc
+	dy = y - yc
+	dX = dx * cos(pa) + dy * sin(pa)
+	dY = (-dx * sin(pa) + dy * cos(pa)) / ar
+	R = sqrt (dX ** 2 + dY ** 2)
+.fi
+
+where dx and dy are the object coordinates relative to the object
+center,  dX and dY are the object coordinates in the transformed
+circular coordinates, and R is the circularly symmetric radius.
+The transformation parameters are the axial ratio \fIar\fR
+defined as the ratio of the minor axis to the major axis,
+and the position angle \fIpa\fR defined counterclockwise from
+the x axis.
+
+The \fIradius\fR parameter defines the size, in pixels, of the model
+object (before seeing for the galaxies) in the output image.  It
+consistently refers to the major axis of the object but its meaning
+does depend on the model.  For the gaussian and moffat profiles it is
+defined as the half-intensity radius.  For the sersic, expdisk, and devauc
+profiles it is defined as the half-flux radius.  For the user specified
+profiles it is the radius of the last profile point.  And for the image
+templates it is the radius of the image along the first or x axis given
+by one-half of the image dimension; i.e. nc/2.
+
+The profiles of the analytic functions extend to infinity so a dynamic
+range, the ratio of the peak intensity to the cutoff intensity, is imposed
+to cutoff the profiles.  The \fIdynrange\fR package parameter applies to
+the stellar and galaxy analytic profiles.  The larger this parameter the
+further the profile extends, particularly for the large index Sersic and De
+Vaucouleurs models.  When modeling large galaxies this has a fairly strong
+affect on the execution time because the overall extent of the images
+becomes rapidly greater.  Only for very high signal-to-noise objects will
+the cutoff be noticeable.  A correction is made to account for lost light
+(light beyond the modeled dynamic range) so that an aperture magnitude
+will give the correct value for an object of the specified total magnitude.
+This can become quite significant for larger index Sersic profiles and
+for the default dynamic range.
+
+The object models are integrated over the size of the image pixels.  This
+is done by subsampling, dividing up a pixel into smaller pieces called
+subpixels.  For the image templates a bilinear surface interpolation
+function is used and integrated analytically over the extent of the
+subpixels.  The user cumulative one dimensional profiles are first
+converted to intensity profiles.  The various intensity profiles are then
+binned into pixel fluxes per subpixel on a grid much finer than the
+subpixel spacing.  Then for any particular radius and object center the
+appropriate subpixel flux can be determined quickly and accurately.
+
+The number of subpixels per image pixel is determined by the package
+parameters \fInxsub\fR, \fInysub\fR, \fInxgsub\fR, and \fInygsub\fR.  The
+first two apply to the stars and the PSF and the latter two apply to the
+galaxies.  Typically the subsampling will be the same in each dimension.
+The galaxies are generally  subsampled less since they will have less
+rapidly changing profiles and are convolved by the PSF.  Also, the stars
+are computed only a few times and then scaled and moved, as described
+below, while each galaxy needs to be computed separately.  Therefore, one
+can afford greater precision in the stars than in the galaxies.
+
+Given an image of several hundred pixels subsampled by a factor of 100
+(10 x 10) this will be a very large number of computations.  A
+shortcut to reduce this number of operations is allow the number
+of subpixels to change as a function of distance from the
+profile center.  Since the profile center is where the intensity
+changes most rapidly with position, the greatest subsampling is needed for
+the pixel nearest the center.  Further from the object center the intensity
+changes more slowly and the number of subpixels may be reduced.
+Thus, the number of subpixels in each dimension in each pixel is
+decreased linearly with distance from the profile center.  For example,
+a pixel which is 3.2 pixels from the profile center will have
+\fInxsub\fR - 3 subpixels in the x dimension.  There is, of course, a
+minimum of one subpixel per pixel or, in other words, no subsampling
+for the outer parts of the objects.  By adjusting the subsampling
+parameters one can set the degree of accuracy desired at the trade off of
+greatly different execution times.
+
+The star shapes are assumed constant over the images and only their
+position and magnitude change.  Thus, rather than compute each desired
+star from the model profile or image template, a normalized star
+template is computed once, using the spatial transformation and
+subsampling operations described above, and simply scaled each time to
+achieve the desired magnitude and added at the requested position.
+However, the apparent star shape does vary depending on where its
+center lies within an image pixel.  To handle this a set of
+normalized star templates is precomputed over a grid of centers
+relative to the center of a pixel.  Then the template with center
+nearest to that requested, relative to a pixel center, is used.  The
+number of such templates is set by the package parameters \fInxc\fR and
+\fInyc\fR where the two axis typically have the same values.  The
+larger the number of centers the more memory and startup time required
+but the better the representation of this sampling effect.  The choice
+also depends on the scale of the stars since the larger the star
+profile compared to a pixel the smaller the subcentering effect is.
+This technique allows generating images with many stars, such as a
+globular cluster or a low galactic latitude field, quite
+efficiently.
+
+Unlike the stars, the galaxies will each have different profiles,
+ellipticities, and position angles and so templates cannot be used (except
+for special test cases as mentioned later).  Another difference is that the
+galaxy models need to be convolved by the PSF; i.e. the shapes are defined
+prior to seeing.  The PSF convolution must also be subsampled and the
+convolution operation requires as many operations as the number of pixels
+in the PSF for each galaxy subpixel.  Thus, computing seeing convolved,
+well subsampled, large galaxy images is the most demanding task of all,
+requiring all the shortcuts described above (larger and variable
+subsampling and the subpixel flux approximation) as well as further ones.
+
+The PSF used for convolving galaxies is truncated at a lower dynamic
+range than the stars according to the package parameter
+\fIpsfrange\fR.  This reduces the number of elements in the convolution
+dramatically at the expense of losing only a small amount of the flux
+in the wings.  Like the stars, the PSF is precomputed on a grid of
+pixel subcenters and the appropriate PSF template is used for each
+galaxy subpixel convolution.  Unlike the stars, the truncated PSF is
+normalized to unit flux in order to conserve the total flux in the
+galaxies.  For the extended galaxies this approximation has only a very
+small effect.  As with the other approximations one may increase the
+dynamic range of the PSF at the expense of an increase in execution
+time.
+
+There is an exception to using the truncated PSF.  If the size of the
+galaxy because very small, 0.01 pixel, then a stellar image is substituted.
+
+
+OBJECT FILES
+
+The object files contain lines defining stars and galaxies.  Stars
+are defined by three numbers and galaxies by seven or eight as
+represented symbolically below.
+
+.nf
+           stars:  xc yc magnitude
+        galaxies:  xc yc magnitude model radius ar pa <save>
+.fi
+
+.ls 6 xc, yc:
+Object center coordinates.  These coordinates are transformed to image
+coordinates as follows.
+
+.nf
+	xc in image = xoffset + xc / distance
+	yc in image = yoffset + yc / distance
+.fi
+
+where \fIxoffset\fR and \fIyoffset\fR are the task offset parameters.
+Objects whose image centers fall outside the image dimensions are ignored.
+.le
+.ls magnitude:
+Object magnitude.  This is converted to instrumental fluxes as follows.
+
+.nf
+	flux = exptime/distance**2 * 10**(-0.4*(magnitude-magzero))
+.fi
+
+where \fIexptime\fR, \fIdistance\fR, and \fImagzero\fR are task parameters.
+For the analytic star and galaxy models a correction
+is made for lost light due to the finite extent of the image in the
+sense that the flux added to the image will never quite be that
+requested.
+.le
+.ls model:
+The types of galaxy models are as follows:
+.ls 4 sersic<n>
+A Sersic model of index n.  The index may real but the value will be rounded
+to the nearest multiple of 0.5 or, equivalently, two times the index value will
+be rounded to an integer.  The index must be between 0.5 and 10.  The Sersic
+model defined as
+
+.nf
+	I = exp (-b * (R/radius)**1/n)
+.fi
+
+where radius is the major axis scale length corresponding to half of the
+total flux.  The value of b is computed using the formula of Ciotti and
+Bertin (AA v352, p447, 1999);
+
+.nf
+	b = 2n - 1/3 + 4/(405n) + 46 / (25515n^2)
+.fi
+.le
+.ls 4 expdisk
+An exponential disk model defined as
+
+.nf
+	I = exp (-b * R/radius)
+.fi
+
+where radius is the major axis scale length corresponding to half of the total
+flux and b is computed as with the Sersic model for n=1.  In fact, the
+algorithm is identical with that for the Sersic model using n=1.  Note that
+because of this there will be slight differences with the earlier versions.
+.le
+.ls devauc
+A De Vaucouleurs profile defined as
+
+.nf
+	I = exp (-b * (R/radius)**1/4)
+.fi
+
+where radius is the major axis scale length corresponding to half of the total
+flux and b is computed as with the Sersic model for n=4.  In fact, the
+algorithm is identical with that for the Sersic model using n=4.  Note that
+because of this there will be slight differences with the earlier versions.
+.le
+.ls <image>
+If not one of the profiles above an image of the specified name is
+sought.  If found the center of the template image is assumed to be the
+center of the object and the image template is scaled so that the
+radius of the template is given by the major axis scale radius parameter.
+.le
+.ls <profile file>
+If not one of the above a text file giving a cumulative flux profile from
+the center to the edge is sought.  If found the profile defines
+a model galaxy of extent to the last profile point given by
+the major axis scale radius parameter.
+.le
+.le
+.ls 6 radius:
+Major axis scale radius parameter in pixels as defined above for the different
+galaxy models.  The actual image radius is modified as follows.
+
+	radius in image = radius / distance
+.le
+.ls ar:
+Minor to major axis axial ratio.
+.le
+.ls pa:
+Major axis position angle in degrees measured counterclockwise from the X axis.
+.le
+.ls save:
+If a large number of identically shaped galaxies (size, axial ratio,
+and position angle) located at the same subpixel (the same x and y
+fractional part) but with varying magnitudes is desired then by
+putting the word "yes" as the eighth field the model will be saved
+the first time and reused subsequent times.  This speeds up the execution.
+There may certain algorithm testing situations where this might be useful. 
+.le
+.ih
+EXAMPLES
+1. Create a galaxy cluster with a power law distribution of field galaxies
+and stars as background/foreground.
+
+.nf
+    ar> gallist galaxies.dat 100 spatial=hubble lum=schecter egal=.8
+    ar> gallist galaxies.dat 500
+    ar> starlist galaxies.dat 100
+    ar> mkobjects galaxies obj=galaxies.dat gain=3 rdnoise=10 poisson+
+.fi
+
+Making the image takes about 5 minutes (2.5 min cpu) on a SPARCstation 1.
+
+2. Create a uniform artificial starfield of 5000 stars for a 512 square image.
+
+.nf
+    ar> starlist starfield.dat 5000
+    ar> mkobjects starfield obj=starfield.dat gain=2 rdnoise=10 poisson+
+.fi
+
+This example takes about a minute on a SPARCstation 1.
+
+3. Create a globular cluster field of 5000 stars for a 512 square image.
+
+.nf
+    ar> starlist gc.dat 5000 spat=hubble lum=bands
+    ar> mkobjects gc obj=gc.dat gain=2 rdnoise=10 poisson+
+.fi
+
+This example takes about a minute on a SPARCstation 1.
+
+4. Add stars to an existing image for test purposes.
+
+.nf
+    ar> mkobjects starfield obj=STDIN gain=2 pois+ magzero=30
+    100 100 20
+    100 200 21
+    200 100 22
+    200 200 23
+    [EOF]
+.fi
+
+5. Look at the center of the globular cluster with no noise and very
+good seeing.
+
+.nf
+	cl> mkobjects gc1 obj=gc.dat nc=400 nl=400 distance=.5 \
+	>>> xo=-313 yo=-313 radius=.1
+.fi
+
+The offset parameters are used to recenter the cluster from
+(256,256) in the data file to (200,200) in the expanded field.
+This example takes 30 sec (5 sec CPU) on a SPARCstation 1.  To expand
+and contract about a fixed point define the object list to have an
+origin at zero.
+
+.nf
+    ar> starlist gc.dat 5000 spat=hubble lum=bands xmin=-256 xmax=256 \
+    >>> ymin=-256 ymax=256
+    ar> mkobjects gc obj=gc.dat xo=257 yo=257 gain=2 rdnoise=10 poisson+
+    ar> mkobjects gc1 obj=gc.dat xo=257 yo=257 gain=2 \
+    >>> distance=.5 rdnoise=10 poisson+
+.fi
+
+6. Make an image of dev$pix at various distances and orientation.  First we
+must subtract the background.
+
+.nf
+	cl> imarith dev$pix - 38 pix
+	cl> mkobjects pix1 obj=STDIN nc=200 nl=200 back=1000 \
+	>>> magzero=30 rd=10 poi+
+	50 50 15.0 pix 40 1 0
+	150 50 15.6 pix 30 .8 45
+	50 150 16.5 pix 20 .6 90
+	150 150 17.1 pix 15 .4 135
+	[EOF]
+.fi
+
+It would be somewhat more efficient to first block average the
+template since the oversampling in this case is very large.
+.ih
+REVISIONS
+.ls MKOBJECTS V2.11+
+The random number seed can be set from the clock time by using the value
+"INDEF" to yield different random numbers for each execution.
+.le
+.ls MKOBJECTS V2.11
+The default value of "ranbuf" was changed to zero.
+.le
+.ih
+SEE ALSO
+gallist, starlist, mknoise, mkheader
+.endhelp
diff --git a/noao/artdata/doc/mkpattern.hlp b/noao/artdata/doc/mkpattern.hlp
new file mode 100644
index 00000000..61bfae09
--- /dev/null
+++ b/noao/artdata/doc/mkpattern.hlp
@@ -0,0 +1,180 @@
+.help mkpattern Jan90 noao.artdata
+.ih
+NAME
+mkpattern - Make/add patterns in images
+.ih
+USAGE
+mkpattern input
+.ih
+PARAMETERS
+.ls input
+Images to create or modify.  Image sections are allowed to apply a pattern
+to a portion of an image.
+.le
+.ls output = ""
+Output images when modifying input images.  If no output images are
+given then existing images in the input list are modified directly.
+If an output image list is given then it must match in number the
+input list.
+.le
+.ls pattern = "constant"
+Pattern to be used.  The patterns are:
+.ls constant
+Constant value v1.
+.le
+.ls grid
+A grid starting with the first pixel and going in steps of the
+pattern size with value v2.  Pixels between the grid have value v1.
+A minimum grid size of 2 is enforced.
+.le
+.ls checker
+A checkerboard with squares of the pattern size alternating between values v1
+and v2 starting with v1.
+.le
+.ls coordinates
+Each pixel is numbered sequentially starting with 1 with the column
+dimension varying fastest.
+.le
+.ls slope
+A sloped plane starting with value v1 for the first pixel and value v2 for
+the last pixel in one or two dimensions.
+.le
+.ls square
+A checkerboard pattern in which the size of the squares begin with
+the pattern size and grow as the square of the coordinate.
+.le
+.le
+.ls option = "replace"
+Editing option when modifying existing images.  Often this is used
+in conjunction with image sections to modify a part of an image.
+The options are:
+
+.nf
+ replace - Replace the image with the pattern.
+     add - Add the pattern to the image.
+multiply - Multiply the pattern with the image values.
+.fi
+.le
+.ls v1 = 0., v2 = 1.
+Pattern values used as described for each pattern.
+.le
+.ls size = 1
+Pattern size used as described for each pattern.
+.le
+
+WHEN CREATING NEW IMAGES
+.ls title = ""
+Image title to be given to the images.  Maximum of 79 characters.
+.le
+.ls pixtype = "real"
+Pixel datatype of new images; one of ushort, short, integer, real, double,
+or complex.
+.le
+.ls ndim = 2
+Number of dimensions between 0 and 7.
+.le
+.ls ncols = 512, nlines = 512
+Number of columns (first dimension) and lines (second dimension).
+.le
+.ls n3 = 1, n4 = 1, n5 = 1, n6 = 1, n7 = 1
+Number of pixels in 3rd-7th  dimensions
+.le
+.ls header = "artdata$stdheader.dat"
+Image or header keyword data file.  If an image is given then the image header
+is copied.  If a file is given then the FITS format cards are copied.
+This only applies to new images.   The data file consists of lines
+in FITS format with leading whitespace ignored.  A FITS card must begin
+with an uppercase/numeric keyword.  Lines not beginning with a FITS
+keyword such as comments or lower case are ignored.  The user keyword
+output of \fBimheader\fR is an acceptable data file.  See \fBmkheader\fR
+for further information.
+.le
+.ih
+DESCRIPTION
+This task creates or modifies images with a choice of patterns.  New images
+are created with the specified dimensions, datatype, and pattern.
+Existing images may have the pattern replace, add, or multiply the
+pixel values.  Existing images may be modified in place or new images may be
+created and image sections are allowed.
+
+For new images a set of header keywords may be added by specifying an
+image or data file with the \fIheader\fR parameter (see also \fBmkheader\fR).
+If a data file is specified lines beginning with FITS keywords are
+entered in the image header.  Leading whitespace is ignored and any
+lines beginning with words having lowercase and nonvalid FITS keyword
+characters are ignored.
+
+This task is the simplest one for creating empty images to be used for
+mosaicing with \fBimcopy\fR and making patterns for testing display and
+image operators.  The replace option is generally used with image sections
+to place constant values in regions.  The multiply option is useful
+for making masks of the given pattern when the values are 0 and 1.
+
+Though the patterns make sense extending to higher dimensions they
+are only defined in two dimensions.  One dimensional images may be
+thought of as the first line of the two dimensional pattern.  Images
+with dimensions greater than 2 simply repeat the two dimensional
+pattern into the higher dimensions.  The reason for stopping at
+two dimensions is simplicity.
+
+The patterns have the following precise definitions where P(i,j) is the
+pixel value at column i and line j, v1 and v2 are the pattern
+values, size is the pattern size, ncols and nlines are the number of
+columns and lines in the image, int is the integer function, mod is the
+modulus function, and sqrt is the square root function.
+
+.nf
+                k = int ((i-1)/size), l = int ((j-1)/size)
+                ksr = int (sqrt (k)), lsr = int (sqrt (l))
+                slope = (v2-v1) / ((ncols+nlines-2)/size)
+
+    constant:   P(i,j) = v1
+
+        grid:   P(i,j) = v2   when mod(i,size)=1 or mod(j,size)=1
+                P(i,j) = v1   otherwise
+
+ coordinates:   P(i,j) = i + j * ncols
+
+     checker:   P(i,j) = v1   when mod(k,2)=0 and mod(l,2)=0
+                P(i,j) = v2   when mod(k,2)=1 and mod(l,2)=0
+                P(i,j) = v2   when mod(k,2)=0 and mod(l,2)=1
+                P(i,j) = v1   when mod(k,2)=1 and mod(l,2)=1
+
+       slope:   P(i,j) = v1 + slope * (k + l) 
+
+      square:   P(i,j) = v1   when mod(ksr,2)=0 and mod(lsr,2)=0
+                P(i,j) = v2   when mod(ksr,2)=1 and mod(lsr,2)=0
+                P(i,j) = v2   when mod(ksr,2)=0 and mod(lsr,2)=1
+                P(i,j) = v1   when mod(ksr,2)=1 and mod(lsr,2)=1
+.fi
+
+.ih
+EXAMPLES
+1. Create an empty (constant value of zero) three dimensional image.
+
+.nf
+	cl> mkpattern cube ndim=3 nc=100 nl=100 n3=100
+.fi
+
+2. Replace a square region of an image with the value -1000.
+
+.nf
+	cl> mkpat alpha[201:250,1:50] v1=-1000
+.fi
+
+3. Put a grid pattern on an image to create a new image.
+
+.nf
+	cl> mkpat dev$pix out=gridpix pat=grid op=mul v1=1 v2=0
+.fi
+.ih
+REVISIONS
+.ls MKPATTERN V2.11
+Now allows ndim=0 to create dataless header.
+
+Now allows type ushort pixel type.
+.le
+.ih
+SEE ALSO
+imcopy, imreplace
+.endhelp
diff --git a/noao/artdata/doc/starlist.hlp b/noao/artdata/doc/starlist.hlp
new file mode 100644
index 00000000..de3eb8b8
--- /dev/null
+++ b/noao/artdata/doc/starlist.hlp
@@ -0,0 +1,355 @@
+.help starlist Feb90 noao.artdata
+.ih
+TASK
+starlist -- make an artificial star list
+.ih
+USAGE
+starlist starlist nstars
+.ih
+PARAMETERS
+.ls starlist
+The name of the output text file for the x and y coordinates
+and magnitudes of the artificial stars.  Output will be appended to this
+file is it exists.
+.le
+.ls nstars = 5000
+The number of stars in the output star list.
+.le
+.ls interactive = no
+Examine plots and change the parameters of the spatial luminosity
+distributions interactively.
+.le
+
+			SPATIAL DISTRIBUTION
+.ls spatial = "uniform"
+Type of spatial distribution.  The types are:
+.ls uniform
+The stars are uniformly distributed between \fIxmin\fR, \fIxmax\fR, \fIymin\fR,
+and \fIymax\fR.
+.le
+.ls hubble
+The stars are distributed around the center of symmetry \fIxcenter\fR and
+\fIycenter\fR according to a Hubble density law of core radius
+\fIcore_radius\fR and background density \fIbase\fR.
+.le
+.ls file
+The radial density function is contained in the text file \fIsfile\fR.
+.le
+.le
+.ls xmin = 1., xmax = 512., ymin = 1., ymax = 512.
+The range of output coordinates in x and y.
+.le
+.ls xcenter = INDEF, ycenter = INDEF
+The coordinate of the center of symmetry for the "hubble"
+and "file" radial density functions. The default is the
+midpoint of the coordinate limits.
+.le
+.ls core_radius = 30
+The core radius of the Hubble spatial distribution in pixels.
+.le
+.ls base = 0.0
+The background density relative to the central density of the Hubble
+density distribution.
+.le
+.ls sseed = 1
+The initial value supplied to the random number generator used to
+generate the output x and y coordinates.
+If a value of "INDEF" is given then the clock
+time (integer seconds since 1980) is used as the seed yielding
+different random numbers for each execution.
+.le
+
+			MAGNITUDE DISTRIBUTION
+.ls luminosity = "powlaw"
+Type of luminosity distribution.  The types are:
+.ls uniform
+The stars are uniformly distributed between \fIminmag\fR and \fImaxmag\fR.
+.le
+.ls powlaw
+The stars are distributed according to a power law with coefficient
+\fIpower\fR.
+.le
+.ls salpeter
+The stars are distributed with a Salpeter luminosity function between
+\fIminmag\fR and \fImaxmag\fR.
+.le
+.ls bands
+The stars are distributed with a Bahcall and Soneira luminosity function
+between \fIminmag\fR and \fImaxmag\fR.  The function is described
+by the parameters \fIalpha\fR, \fIbeta\fR, \fIdelta\fR and \fImstar\fR
+whose default values give a best fit to the observed main sequence in several
+nearby globular clusters.
+.le
+.ls file
+The luminosity function is contained in the text file \fIlfile\fR.
+.le
+.le
+.ls minmag = -7., maxmag = 0.
+The range of output magnitudes.  The "salpeter" luminosity function
+imposes limits of -4 and 16 and the "bands" luminosity function
+imposes limits of -7 and 17 relative to the zero point given by
+\fImzero\fR.
+.le
+.ls mzero = -4.
+The zero point for converting the output relative magnitudes
+to absolute magnitudes for the Salpeter and Bahcall and Soneira
+luminosity functions.  For example the default values give an
+absolute magnitude range of -3 to +4.
+.le
+.ls power = 0.6
+Coefficient for the power law magnitude distribution.
+The default value of 0.6 is the value for a homogeneous
+and isotropic distribution with no cutoff in distance.
+.le
+.ls alpha = 0.74, beta = 0.04, delta = 0.294, mstar = 1.28
+The parameters of the Bahcall and Soneira luminosity function.
+.le
+.ls lseed = 1
+The initial value supplied to the random number generator used to
+generate the output magnitudes.
+If a value of "INDEF" is given then the clock
+time (integer seconds since 1980) is used as the seed yielding
+different random numbers for each execution.
+.le
+
+			USER FUNCTIONS
+.ls sfile
+The name of the input text file containing the sampled spatial radial
+density
+function, one sample point per line, with the radius and relative probability
+in columns one and two respectively. The sample points need not be
+uniformly spaced or normalized.
+.le
+.ls nssample = 100
+The number of points at which the \fIspatial\fR density function is 
+sampled. If the \fIspatial\fR density function is analytic or approximated
+analytically (the "uniform" and "hubble" options) the function is sampled
+directly. If the function is read from a file  (the "file" option) an
+initial smoothing step is performed before sampling.
+.le
+.ls sorder = 10
+The order of the spline fits used to evaluate the integrated spatial
+density function.
+.le
+.ls lfile
+The name of the input text file containing the sampled luminosity
+function, one sample point per line, with the magnitude and relative probability
+in columns one and two respectively. The sample points need not be
+uniformly spaced or normalized.
+.le
+.ls nlsample = 100
+The number of points at which the luminosity function is sampled. If
+the luminosity function is analytic or approximated analytically (the
+"salpeter" and "bands" options) the function is sampled directly.  If
+it is read from a file  (the "file" option) an initial smoothing step
+is performed before sampling.
+.le
+.ls lorder = 10
+The order of the spline fits used to evaluate the integrated
+\fIluminosity\fR function.
+.le
+
+			INTERACTIVE PARAMETERS
+.ls rbinsize = 10.
+The bin size in pixels of the plotted histogram of the radial density
+distribution.
+.le
+.ls mbinsize = 0.5
+The bin size in magnitudes of the plotted histogram of the luminosity function.
+.le
+.ls graphics = stdgraph
+The default graphics device.
+.le
+.ls cursor = ""
+The graphics cursor.
+.le
+.ih
+DESCRIPTION
+\fBStarlist\fR generates a list of x and y coordinates and magnitudes
+for a sample of \fInstars\fR stars based on a user selected spatial
+density function \fIspatial\fR  and luminosity function
+\fIluminosity\fR and writes (appends) the results to the text file
+\fIstarlist\fR. If the \fIinteractive\fR parameter is "yes" the user
+can interactively examine plots of the spatial density function,
+the radial density function, and the luminosity function, and alter the
+parameters of the task until a satisfactory artificial field is
+generated.
+
+The spatial density function generates x and y values around a center
+of symmetry defined by \fIxcenter\fR and \fIycenter\fR within the x and
+y limits \fIxmin\fR, \fIxmax\fR, \fIymin\fR and \fIymax\fR according to
+the spatial density function specified by \fIspatial\fR.  The three
+supported spatial density functions are listed below where R is the
+radial distance in pixels, P is the relative spatial density, C is a
+constant and f is the best fitting cubic spline function to the spatial
+density function R(user), P(user) supplied by the user in the text file
+\fIsfile\fR.
+
+.nf
+    uniform:  P = C
+    hubble:   P = 1.0 / (1 + R / core_radius) ** 2 + base
+    file:     P = f (R(user), P(user))
+.fi
+
+The Hubble and user file spatial density function are sampled at
+\fInssample\fR equally spaced points, and integrated to give the
+spatial density probability function at each sampled point. The
+integrated probability function is normalized and approximated by a
+cubic spline of order \fIsorder\fR.  The x and y coordinates are
+computed by randomly sampling the integrated probability function until
+\fInstars\fR stars which satisfy the x and y coordinate limits
+\fIxmin\fR, \fIxmax\fR, \fIymin\fR and \fIymax\fR are generated.
+
+The luminosity function generates relative magnitude values between
+\fIminmag\fR and \fImaxmag\fR according to the luminosity function
+specified by \fIluminosity\fR.  The four supported luminosity functions
+are defined below where M is the magnitude, P is the relative luminosity
+function, C is a constant and f is the best fitting cubic spline
+function to the luminosity function M(user), P(user) supplied by the
+in the text file \fIlfile\fR.
+
+.nf
+  uniform:  P = C
+
+  powlaw:   P = C * 10. ** (power * M)
+
+  salpeter: P = C * 10. ** (-3.158 + 1.551e-1*dM - 5.194e-3*dM**2)
+
+            dM = M - mzero
+
+                             C * 10. ** (beta * dM)
+  bands:   P =  --------------------------------------------------
+               (1. + 10. ** ((beta-alpha)*delta*dM))) ** 1. /delta
+
+           dM = M - mstar - mzero
+
+  file:    P = f (M(user), P(user))
+.fi
+
+The Salpeter and "bands" functions are defined in terms of absolute
+magnitudes so the parameter \fImzero\fR is used to convert from
+relative magnitudes.  Equivalently, one could use absolute magnitudes
+for the magnitude limits while setting the zero point to 0.
+
+The luminosity function is sampled at \fInlsample\fR equally spaced
+points, and integrated to give the luminosity probability function at
+each sampled point. The probablity function is normalized and
+approximated by a cubic spline of order \fIlorder\fR. The magnitudes
+are computed by randomly sampling the integrated probability function
+until \fInstars\fR objects which satisfy the magnitude limits
+\fIminmag\fR and \fImaxmag\fR are generated.  The Salpeter luminosity
+is a best fit function to the data of McCuskey (McCuskey, 1966, Vistas
+Astr. 7, 141). The Bahcall and Soneira function and the default values
+of the parameters are discussed by Bahcall and Soneira (Ap.J.  Supp. 44, 73).
+.ih
+CURSORS
+The following interactive keystroke commands are available from within the
+STARLIST task.
+
+.nf
+	Starlist Keystroke Commands
+
+?	Print options
+f	Fit  one or more of the following
+	    Spatial density function (SDF)
+	    Luminosity functions (LF)
+x	Plot the x-y spatial density function
+r	Plot the histogram of the radial density function
+m	Plot the histogram of the luminosity function
+:	Colon escape commands (see below)
+q	Exit program
+.fi
+
+The following parameters can be shown or set from within the STARLIST task.
+
+
+.nf
+		Starlist Colon Commands
+
+:show			Show starlist parameters
+:nstars     [value]	Number of stars
+
+:spatial    [string]	Spatial density function (SDF)
+			(uniform|hubble|file) 
+:xmin       [value]	Minimum X value
+:xmax       [value]	Maximum X value
+:ymin       [value]	Minimum Y value
+:ymax       [value]	Maximum Y value
+:xcenter    [value]	X center for SDF
+:ycenter    [value]	Y center for SDF
+:core       [value]	Core radius for Hubble density function
+:base       [value]	Background density for Hubble density function
+
+:luminosity [string]	Luminosity function (LF)
+			(uniform|powlaw|salpeter|bands|file)
+:minmag     [value]	Minimum magnitude
+:maxmag     [value]	Maximum magnitude
+:mzero	    [value]	Magnitude zero-point for salpeter and bands LF
+:power	    [value]	Exponent for powlaw LF
+:alpha      [value]	Alpha parameter for bands LF
+:beta       [value]	Beta parameter for bands LF
+:delta      [value]	Delta parameter for bands LF
+:mstar      [value]	Mstar parameter for bands LF
+
+:sfile	    [string]    File containing the user SDF
+:nssample   [value]	Number of SDF sample points
+:sorder	    [value]	Order of spline fit to integrated SDF
+:lfile	    [string]    File containing the user LF
+:nlsample   [value]	Number of LF sample points 
+:lorder	    [value]	Order of spline fit to the integrated LF
+
+:rbinsize   [value]	Resolution of radial profile histogram (pixels)
+:mbinsize   [value]	Resolution of magnitude histogram (mag)
+.fi
+
+.ih
+EXAMPLES
+1. Create a uniform artificial starfield of 5000 stars for a 512 square image.
+
+.nf
+    ar> starlist starfield.dat 5000
+    ar> mkobjects starfield obj=starfield.dat gain=2 rdnoise=10 poisson+
+.fi
+
+This example takes about a minute on a SPARCstation 1.
+
+2. Create a globular cluster field of 5000 stars for a 512 square image.
+
+.nf
+    ar> starlist gc.dat 5000 spat=hubble lum=bands
+    ar> mkobjects starfield obj=gc.dat gain=2 rdnoise=10 poisson+
+.fi
+
+This example takes about a minute on a SPARCstation 1.
+
+3. Examine the distributions for a Hubble spatial distribution
+and Salpeter magnitude distribution using 1000 stars without
+creating a data file.
+
+.nf
+    ar> starlist dev$null 1000 inter+ spat=hubble lum=salpeter
+	    ... an x-y plot will appear on the screen
+	    ... type r to examine the radial density function
+	    ... type m to examine the luminosity function
+	    ... type = to make a copy of any of the plots
+	    ... type q to quit
+.fi
+.ih
+REVISIONS
+.ls STARLIST V2.11+
+The random number seeds can be set from the clock time by using the value
+"INDEF" to yield different random numbers for each execution.
+.le
+.ih
+BUGS
+The spline approximation to the spatial density and luminosity
+probability functions can  cause wiggles in the output spatial density
+and luminosity functions. Users can examine the results interactively
+and experiment with the spline order and number of sample points if
+they are not satisfied with the results of STARLIST. The default setup
+of 10 sample points per spline piece is generally satisfactory for the
+spatial density and luminosity functions supplied here.
+.ih
+SEE ALSO
+gallist mkobjects
+.endhelp
diff --git a/noao/artdata/doc/version1.1 b/noao/artdata/doc/version1.1
new file mode 100644
index 00000000..b16c4df1
--- /dev/null
+++ b/noao/artdata/doc/version1.1
@@ -0,0 +1,49 @@
+		Changes to the ARTDATA Package: Version 1.1
+
+Since the first release of the ARTDATA package with V2.9 there have been
+several bugs found and documented, as is natural, and some new features
+and tasks have been added.  This note summarizes these changes.
+
+The bugs (log numbers 127, 130, 131, and 135) are all in the main task
+MKOBJECTS.  Two of these concern the size of the PSF as specified by
+the "radius" parameter.  The radius parameter for the gaussian PSF
+is interpreted as a FWHM instead of the intended HWHM (radius at half
+maximum).  This leads to very narrow stars until one compensates for the
+error.  The size of the moffat PSF are off by a small factor (~10%)
+depending on the beta parameter.  Another bug causes objects near the edge
+to be off by 1 pixel from the specified coordinates.  The last bug
+prevents use of a user supplied profile function for the PSF.  A
+distributed patch for the above problems is planned.  Thanks go to
+Lindsey Davis and Steven Perry for the careful studies which identified
+the first three bugs.
+
+The artificial one dimensional spectra produced by MK1DSPEC have new
+parameters allowing addition of a velocity shift specified either as
+a velocity or redshift.  This is a useful addition for tests of radial
+velocity programs.
+
+The initial version of the package provided minimal headers; mostly
+just the wavelengths for 1D spectra and the dispersion axis for 2D
+spectra.  People then need to HEDIT keywords to make the images
+work well with some packages.  Four things have been done to
+make useful headers for the artificial data.  First, a detailed
+log of parameters used to generate the images is included as COMMENT
+cards.  This includes task parameters and information from data files
+such as those in the star and galaxy list files produced by
+STARLIST and GALLIST.  Second, the gain, read noise, exposure
+time, and dispersion correction flag have been added to the header
+when appropriate.  Third, each task has a new parameter for specifying
+a header keyword data file containing a list of keywords and values
+to be automatically added.  A default file is supplied but any set
+for a particular type of data may be substituted.  Finally, a new
+task called MKHEADER has been added which applies header keyword
+data files to images.
+
+A new task called MKEXAMPLES has been added.  Given the name of an
+example from a menu an image is created.  This task is intended to
+provide examples for the ARTDATA package as well as test and
+demonstration images for the various packages.  The initial menu
+includes long slit and multifiber spectra, a globular cluster,
+a star field, a galaxy field, and a galaxy cluster.  Additional
+examples will be added as demonstrations and test procedures are
+developed for various packages.