;+
; NAME:
;   GAUSS1
;
; AUTHOR:
;   Craig B. Markwardt, NASA/GSFC Code 662, Greenbelt, MD 20770
;   craigm@lheamail.gsfc.nasa.gov
;
; PURPOSE:
;   Compute Gaussian curve given the mean, sigma and area.
;
; MAJOR TOPICS:
;   Curve and Surface Fitting
;
; CALLING SEQUENCE:
;   YVALS = GAUSS1(XVALS, [MEAN, SIGMA, AREA], SKEW=skew)
;
; DESCRIPTION:
;
;  This routine computes the values of a Gaussian function whose
;  X-values, mean, sigma, and total area are given.  It is meant to be
;  a demonstration for curve-fitting.
;
;  XVALS can be an array of X-values, in which case the returned
;  Y-values are an array as well.  The second parameter to GAUSS1
;  should be an array containing the MEAN, SIGMA, and total AREA, in
;  that order.
;
; INPUTS:
;   X - Array of X-values.
;
;   [MEAN, SIGMA, AREA] - the mean, sigma and total area of the 
;                         desired Gaussian curve.
;
; INPUT KEYWORD PARAMETERS:
;
;   SKEW - You may specify a skew value.  Default is no skew.
;
; RETURNS:
;
;   Returns the array of Y-values.
;
; EXAMPLE:
;
;   p = [2.2D, 1.4D, 3000.D]
;   x = dindgen(200)*0.1 - 10.
;   y = gauss1(x, p)
;
;   Computes the values of the Gaussian at equispaced intervals
;   (spacing is 0.1).  The gaussian has a mean of 2.2, standard
;   deviation of 1.4, and total area of 3000.
;
; REFERENCES:
;
; MODIFICATION HISTORY:
;   Written, Jul 1998, CM
;
;-

function gauss1, x, p, skew=skew, _EXTRA=extra

  z = (x-p(0))/p(1)
  f = z * 0

  ;; IDL will report "underflow" errors unless the domain is
  ;; restricted.  Here I restrict to 10-sigma
  wh = where(abs(z) LT 10., ct)

  ;; This is a gaussian whose total area is unity
  if ct GT 0 then $
    f(wh) = exp(-0.5D * z(wh)^2)/(sqrt(2.D * !dpi)*p(1))

  ;; Add skew if necessary
  if n_elements(skew) GT 0 then $
    f = (1.D + skew * z)*f

  ;; Apply normalization if needed
  if n_elements(p) GE 3 then f = p(2)*f
  return, f
end