c
c##############################################################################
c
c
c###################################################################
c
      subroutine fnagex2(pt, pcf, kncf, pfit, kermsg)
c
c     this is an iaea subroutine to calculate the electron impact
c     excitation cross sections as a function of the electron
c     energy in eV.
c
c     Cross sections are calculated from the collision strengths:
c	xs (pi*a0**2) = cs / wi / Ee (Ryd) where
c		Ee energy of the impact electron in Rydberg
c		wi statical weight of the initial state
c
c     pt = electron energy in ev
c     pcf is the coefficient data array, where
c     pcf(1)  = itype, index for type of fit, either 1 or 2
c     pcf(2)  = excitation energy , v (also referred to as delta e)
c     pcf(3)  = lower limit of fitting of reduced energy x
c     pcf(4)  = upper limit of fitting of reduced energy x
c     pcf(5)  = statistical weight of initial state (2s+1)*(2l+1)
c     pcf(6)  = parameter A
c     pcf(7)  = parameter B
c     pcf(8)  = parameter C
c     pcf(9)  = parameter D
c     pcf(10) = parameter E
c	type 1:	cs = A + B/X + C/(X*X) + D/(X*X*X) + E*log(X)
c
c if ftype = 2, the pcf array element 10 can be followed by
c     pcf(11)  = parameter F
c	type 2:	cs = A/(X*X) + B*exp(-F*X) + C*exp(-2*F*X) + D*exp(-3*F*X) + E*exp(-4*F*X)
c
c This can be followed by three more parameters for
c the region which contains resonances represented by a linear term
c 	resonnance region: cs = PX + Q where 1 < X < X1
c type 1:
c     pcf(11)  = parameter P
c     pcf(12)  = parameter Q
c     pcf(13)  = parameter X1, the upper limit of the range over which
c                the collision strength is represented by a linear
c                approximation.
c type 2:
c     pcf(12)  = parameter P
c     pcf(13)  = parameter Q
c     pcf(14)  = parameter X1, the upper limit of the range over which
c                the collision strength is represented by a linear
c                approximation.
c
c     kermsg = blank if no errors
c     pfit = cross section in cm2
c
c     written by D. Humbert , iaea atomic and molecular data unit
c     fit function from ADNDT,33,149,(1985) y. itikawa et al, nagoya
c     22 August 2006
c
c------------------------------------------------------------------------
c
      double precision pt, pcf, pfit
      dimension pcf(14)
      character*(*) kermsg
      data s/0.5/
c
      kermsg = ' '
      itype = pcf(1)
      pia02 = .87973554e-16
      rydberg = 13.6056981
      wi = pcf(5)
      vif = pcf(2)
      A = pcf(6)
      B = pcf(7)
      C = pcf(8)
      D = pcf(9)
      E = pcf(10)
c
      X = pt/vif
c
      if  (itype .eq. 2) then
          F = pcf(11)
	  ioffset = 1
      else
          ioffset = 0
      endif
      if (kncf .ge. 12) then
          P  = pcf(11 + ioffset)
          Q  = pcf(12 + ioffset)
          X1 = pcf(13 + ioffset)
          itype = 3
      else
	  X1 = 0.
      endif
c
c---    rate coefficient with resonances: itype = 3
      if (itype .eq. 3) then
         cs = P*X + Q
c---    rate coefficient without resonances
      else if (itype .eq. 1) then
        cs = A + B/X + C/(X*X) + D/(X*X*X) + E*log(X)
	else if (itype .eq. 2) then
	  cs = A/(X*X) + B*exp(-F*X) + C*exp(-2*F*X) +
     1	       D*exp(-3*F*X) + E*exp(-4*F*X)
          else
           kermsg = ' invalid integer for fit type in fnagex (1 or 2)'
           return
      endif
c
c     collision strengh to cross section
      pfit = (cs * pia02 * rydberg) / (wi * pt)
c
      if (pfit .le. 0.0d0) then
            kermsg =
     1  ' error reaction rate is negative check data and temp. range '
          return
        endif
c
      end

def fnagex2(pt, pcf, kncf):
    """
    This function calculates the electron-impact excitation cross sections
     as a function of the  electron energy in eV.

    pt: electron energy in eV
    pcf: parameter data array
        pcf[0]: type of the fit,  either 1 or 2
        pcf[1]: excitation energy
        pcf[2]: lower limit of fitting of reduced energy x
        pcf[3]: upper limit of fitting of reduced energy x
        pcf[4]: statistical weight of initial state (2s+1)*(2l+1)
        pcf[5:14]: fitting parameters, depending upon the type of fit
    """
    itype = pcf[0]
    pia02 = 0.87973554e-16
    rydberg = 13.6056981
    wi = pcf[4]
    vif = pcf[1]
    A = pcf[5]
    B = pcf[6]
    C = pcf[7]
    D = pcf[8]
    E = pcf[9]
    
    X = pt / vif

    if itype == 2:
        F = pcf[10]
        ioffset = 1
    else:
        ioffset = 0

    if kncf >= 12:
        P = pcf[11 + ioffset]
        Q = pcf[12 + ioffset]
        X1 = pcf[13 + ioffset]
        itype = 3
    else:
        X1 = 0.0

    if itype == 3:
        cs = P * X + Q
    elif itype == 1:
        cs = A + B / X + C / (X * X) + D / (X * X * X) + E * np.log(X)
    elif itype == 2:
        cs = A / (X * X) + B * np.exp(-F * X) + C * np.exp(-2 * F * X) + \
             D * np.exp(-3 * F * X) + E * np.exp(-4 * F * X)
    else:
        raise ValueError("Invalid integer for fit type in fnagex (1 or 2)")
    
    pfit = (cs * pia02 * rydberg) / (wi * pt)

    if pfit <= 0.0:
        raise ValueError("Error: Reaction rate is negative. Check data and temperature range.")
    
    return pfit