
      SUBROUTINE gcvmkl(ear, ne, param, ifl, photar, photer)

      INTEGER ne, ifl
      REAL ear(0:ne), param(19), photar(ne), photer(ne)

c
c XSPEC model subroutine to calculate a Gaussian Continous Emission
c Measure distribution.
c   Q(T) = Norm*exp((T-Tmean)^2/2/sigma^2)
c
c This program calls 'vmekal' thus allowing one to vary the elemental
c abundances.
c 
c Parameters:
c    param(1) = mean temperature (keV)
c    param(2) = sigma on temperature
c    param(3) = nH (cm^-3)  Fixed at 1 for most applications
c    param(4) = He abundance
c    param(5) = C   "
c    param(6) = N   "
c    param(7) = O   "
c    param(8) = Ne  "
c    param(9) = Na  "
c    param(10)= Mg  "
c    param(11)= Al  "
c    param(12)= Si  "
c    param(13)= S   "
c    param(14)= Ar  "
c    param(15)= Ca  "
c    param(16)= Fe  "
c    param(17)= Ni  " 
c    param(18) = redshift used in the Mewe-Kaastra plasma model
c    param(19) = switch(0=calculate MEKAL model, 1=interpolate MEKAL model)


      INCLUDE 'xspec.h'

      INTEGER itarr, idem


      REAL sigma, tmean, temp, alow, ahi, winv, tmin, tdelta

      INTEGER switch, i, nt, ntsave, ierr

      CHARACTER contxt*127

      SAVE itarr, idem, ntsave

      REAL erf
      EXTERNAL erf

      DATA itarr, idem, ntsave/0, 0, -1/

      ierr = 0

c Initialize values:

      tmean = param(1)
      sigma = MAX(0., param(2))

c Trap the special case of zero sigma

      IF ( sigma .EQ. 0. ) THEN

         nt = 1
         MEMR(itarr) = tmean
         MEMR(idem) = 1.0

      ELSE

         winv = 1./sigma/sqrt(2.)

c We will integrate over the distribution in 0.2 sigma steps and go out to
c the 3 sigma point in both directions unless that takes the temperature
c below zero.

         tdelta = sigma/5.
  
         tmin = tmean - 3 * sigma
         DO WHILE ( tmin .LE. tdelta/2. )
            tmin = tmin + tdelta
         ENDDO

         nt = 1 + 30 + MIN(30,INT((tmean-tmin)/(tdelta)))

         IF ( nt .NE. ntsave ) THEN
            CALL udmget(nt, 6, itarr, ierr)
            contxt = 'GCVMKL: Failed to get memory for tarr'
            IF ( ierr .NE. 0 ) GOTO 999
            CALL udmget(nt, 6, idem, ierr)
            contxt = 'GCVMKL: Failed to get memory for dem'
            IF ( ierr .NE. 0 ) GOTO 999
            ntsave = nt
         ENDIF

         DO i = 1, nt

            temp = tmin + (i-1)*tdelta
 
            IF ( i .EQ. 1 ) THEN
               alow = 0.5*erf(winv*(temp-tmean-tdelta/2))
               ahi  = 0.5*erf(winv*(temp-tmean+tdelta/2))
            ELSE
               alow = ahi
               ahi  = 0.5*erf(winv*(temp-tmean+tdelta/2))
            ENDIF

            MEMR(itarr+i-1) = temp
            MEMR(idem+i-1)  = ahi - alow

         ENDDO

      ENDIF

      switch = NINT(param(19))

      CALL sumdem(2, switch, ear, ne, param(4), param(5), param(18),
     &            nt, MEMR(itarr), MEMR(idem), ifl, photar, ierr)
      contxt = 'GCVMKL: failure in SUMDEM'

 999  CONTINUE
      IF ( ierr .NE. 0 ) THEN
         CALL xwrite(contxt, 10)
         WRITE(contxt, '(a,i5)') 'Error number = ', ierr
         CALL xwrite(contxt, 10)
      ENDIF

      return
      end 












