c*** File tdclib/specsmooth.f
c*** April 27, 1993
c*** By Doug Mink, Center for Astrophysics

c--- Smooth data from a TDC reduced spectrum file or files

        Subroutine SPECSMOOTH (PARMS,DATASET,NPIX,NEWDATA)

	Character PARMS*(*)
c				String with ASCII parameter value assignments
	Integer*4 DATASET
c				Data set to plot (default should be 1)
	Integer*4 NPIX
c				Number of pixels in spectrum
	Real*4 NEWDATA(*)
c				Smoothed data (returned)

	Integer*4 NDSAMP
c				No. of spectrum points averaged before plotting
	Character*20 NDMODE
c				Processing (sample | average | smooth)
	Character*10 NDTYPE
c				Type of smoothing (gaussian | hanning)
	Integer*4 NDSAMP2
c				Half of NDSAMP
	Real*8 NDHALF
c				Half width of Gaussian in pixels
	Real*8 WT(100)
c				Gaussian weighting function per pixel
	Real*8 WTPIX
c				Gaussian weighting function total

	Logical HANNING
	Real*8 TWOPI,SIG2,RAD2,SUMDATA
	Real*4 SPJ,DATA0,DATA1,DATA2
	Integer*4 IPIXN
	Integer*4 ISAMP,I,IPIX,JWT,JPIX,JPIX0,JPIX1

	Include 'wpconv.cmb'
	Include 'plotdata.cmb'

	NDMODE = 'sample'
	Call PGETC (PARMS,'NDMODE',NDMODE)

	Do 100 IPIX = 1, NPIX
100	    NEWDATA(IPIX) = SPECDATA(IPIX,DATASET)

	If ((NDMODE(1:3).ne.'smo').and.(NDMODE(1:3).ne.'SMO')) Return

	Call PGETI4 (PARMS,'NDSAMP',NDSAMP)
	NDSAMP2 = NDSAMP / 2
	NDTYPE = 'hanning'
	Call PGETC (PARMS,'NDTYPE',NDTYPE)
	If ((NDTYPE(1:1) .eq. 'H') .or. (NDTYPE(1:1) .eq. 'h')) Then
	    HANNING = .true.
	    TWOPI = 6.2831853
	Elseif ((NDTYPE(1:1) .eq. 'N') .or. (NDTYPE(1:1) .eq. 'n')) Then
	    Do 110 I =1,NDSAMP
		DATA0 = NEWDATA(1)
		DATA1 = NEWDATA(1)
		Do 105 IPIX = 1, NPIX
		    IPIXN = IPIX + 1
		    If (IPIXN .gt. NPIX) IPIXN = NPIX
		    DATA2 = NEWDATA(IPIXN)
		    NEWDATA(IPIX) = (DATA0 + DATA1 + DATA1 + DATA2) / 4.0
		    DATA0 = DATA1
105		    DATA1 = DATA2
110		Continue
	    Return
	Else
	    HANNING = .false.
	    Call PGETR8 (PARMS,'NDHALF',NDHALF)
	    SIG2 = 2.d0 * NDHALF * NDHALF
	Endif

c  Set up weighting functions for Gaussian and Hanning smoothing

	NDSAMP = 1 + (NDSAMP2 * 2)
	Do 190 ISAMP = -NDSAMP2, NDSAMP2
	    I = ISAMP + NDSAMP2 + 1
	    RAD2 = Dble (ISAMP) * Dble (ISAMP)
	    If (HANNING) Then
		WT(I) = 0.5 * (1.0 - Dcos (TWOPI*Dble(I) / Dble(NDSAMP-1)))
	    Else
		WT(I) = Exp (-RAD2 / SIG2)
	    Endif
190	    Continue

	Do 300 IPIX = 1, NPIX
	    JPIX0 = IPIX - NDSAMP2
	    JPIX1 = IPIX + NDSAMP2
	    SUMDATA = 0.d0
	    WTPIX = 0.d0
	    JWT =0
	    Do 200 JPIX = JPIX0, JPIX1
		JWT = JWT + 1
		If ((JPIX .ge. 1) .and. (JPIX .le.  NPIX)) Then
		    SUMDATA = SUMDATA + WT(JWT)*SPECDATA(JPIX,DATASET)
		    WTPIX = WTPIX + WT(JWT)
		    Endif
200		Continue
	    If (WTPIX .gt. 0) Then
		SPJ = Sngl (SUMDATA / WTPIX)
	    Else
		SPJ = 0.
	    Endif
300	    NEWDATA(IPIX) = SPJ

	Return
	End

c Feb  6 1991	Add nova-style smoothing
c May 23 1991	Declare undeclared vaariables

c Apr 27 1993	Change double quotes to single quotes
