; NAME:
;      SPECLINES
;
; PURPOSE: 
;      Create a 1D array with speclines
;
; CALLING SEQUENCE:
;      Result = SPECLINES(Header, Maxval, Numl, Ty, Inspec=inspec,
;      fromoh=fromoh)
;
; INPUTS:
;      Header  = The header for the spectrum created
;      Baselvl = The base level for the spectrum
;      Maxval  = The maxium value the lines should have
;      Numl    = The number of lines to be createed
;      Ty      = Size in y direction
;
; KEYWORDS:
;      Inspec  = Contains the user inputed spectrum
;      FromOh  = Set to take lines out of oh line database
;      Bbcurve = Set to use a black body curve as the base
;      Numabsl = Contains the number of absorbtion lines in spectrum
;      Specoff = Contains the object spectrum's offset
;      NoLines = Set to not add any lines
;      OHLoc =   The location of the oh line database
;
; OUTPUT:
;      A 1D array with speclines
;      
; PROCEDURE:
;      If needed create a bbcurve, or else just use a flat baseline.
;      Then if an input spectrum is specified expand it to length ty.
;      Else, if the lines are to be created from the ohdata base, call
;      ohlines, and scale the result.  If random lines are to be used,
;      call randlines and scale.  If numabsl is specified also add
;      absorbtion lines.  Finally convol the result with a gaussing to
;      mimick spectral resolution. 
;
; FUNCTIONS CALLED:
;      sxpar, planckfunc, ohlines, randlines, cnvlgauss
;
; MODIFICATION HISTORY:
;      created Augest 1 2003 by John Dermody

function SPECLINES, header, baselvl, maxval, numl, ty, inspec=inspec, $
                    fromoh=fromoh, bbcurve=bbcurve, numabsl=numabsl, $
                    specoff=specoff, nolines=nolines, ohloc=ohloc

specres = sxpar(header, 'SPECRES')
swavel  = sxpar(header, 'SWAVEL')
ewavel  = sxpar(header, 'EWAVEL')
smile   = sxpar(header, 'SMILE')
ny =      sxpar(header, 'NAXIS2')
by = float(ty - ny) / 2
wavl = (findgen(ty) - by + 1) / ny * $       ; = Wavelength array for CCD
       (ewavel - swavel) + swavel > 0 

; set base level
if (keyword_set(bbcurve)) then begin
  temp = sxpar(header, 'TEMP')
  bbvector = planckfunc(wavl, temp)
  spec = baselvl * bbvector / max(bbvector)
  sxaddhist, 'Object Planck function computed', header
endif else begin
  spec = fltarr(ty) + baselvl
  sxaddhist, 'Set background level', header
endelse

; add spectrum lines
if not keyword_set(nolines) then begin
  if (keyword_defined(inspec)) then begin
    mspec = inspec
    ny = (size(mspec))[1]
    if (ny ne ty) then mspec = [fltarr(by + specoff), mspec]
    ny = (size(mspec))[1]
    if (ny ne ty) then mspec = [mspec, fltarr(ty - ny)]
    spec = spec + mspec / max(mspec) * maxval
    sxaddhist, 'Spectrum from user input', header
  endif else if (keyword_set(fromoh)) then begin
    ohlist = ohlines(wavl, ohloc=ohloc)
    ohspec = addlines(ohlist, ty)
    if (where(ohspec gt 0) ne [-1]) then begin
      spec = spec + ohspec / max(ohspec) * maxval       ; scale
    endif
    sxaddhist, 'Sky spectrum from OH Lines database', header
  endif else begin
    emilist = randlines(ty, numl*100)
    spec = spec + maxval * addlines(emilist, ty)
    if keyword_defined(numabsl) then begin
      abslist = randlines(ty, numabsl*100)
      spec = spec * (1. - addlines(abslist, ty))
    endif
    sxaddhist, 'Randomly generated spectrum', header
  endelse
endif

; smooth from spectral resolution
spec = cnvlgauss(spec, fwhm=specres)

return, spec
end